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In this review article I attempt to summarise past and present-ongoing-work on the problem of 
the inspiral of a small body in the gravitational field of a much more massive Kerr black hole. Such 
extreme mass ratio systems, expected to occur in galactic nuclei, will constitute prime sources of 
gravitational radiation for the future LISA gravitational radiation detector. The article's main goal is 
to provide a survey of basic celestial mechanics in Kerr spacetime and calculations of gravitational 
waveforms and backreaction on the small body's orbital motion, based on the traditional 'flux- 
balance' method and the Teukolsky black hole perturbation formalism. 
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! I. INTRODUCTION 

O . 

■ The launch and subsequent operational deployment of the LISA space-based gravitational wave detector 0,0 is 
an event eagerly awaited by the relativistic astrophysics community. Unlike the various ground-based detectors (both 

^ '. existent and scheduled for the near future) which will be capable of 'scoring' only marginal detections of gravitational 

■ radiation sources, LISA is designed as a prime instrument with the expectation to provide us with unprecedented 
data on strong-gravity spacetime dynamics. 

Mounting astronomical observations point to the existence of massive 'dark objects' in the majority of galactic 

^ ^ , nuclei (see Ref. (3| for an excellent review). Moreover, accretion disc models crucially depend on the existence of an 

J> ■ event- horizon in order to explain the 'dimness' of galactic nuclei such as our own Milky Way 0. It is universally 
' accepted that in reality these objects are supermassive Kerr black holes. 

. Such black holes would strongly interact with their surrounding 'cusp' stellar population, and occasionally capture 
' small-mass /i ^ 1 — 10 Mq compact bodies (black holes, neutron stars and white dwarfs) that happened to scatter 
towards their vicinity Following the capture, the small body will suffer loss of energy and angular momentum by 
emitting gravitational radiation and slowly inspiral until its final plunge into the black hole. The dynamics of such 
\ extreme-mass-ratio-inspiral (EMRI) is (almost) purely gravitational, accurately modelled as a 'test-particle' moving 
in a Kerr spacetime. Due to the extreme mass ratio (the massive black hole should have a mass M ~ IQ^Mq in 
order for the emitted signal to fall in LISA's bandwidth), the evolution of the system is adiabatic; for many orbital 
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I . periods, the small body moves in an almost geodesic trajectory of the background spacetime. Radiation reaction 
^jPj" effects become noticeable at much longer timescales (~ M^//i) and the precise calculation of this motion is one of the 
. . , main current research objectives. 

^ ■ EMRIs, being relatively simple systems, will perhaps prove to be among the 'cleanest' sources for LISA, making 
them unique probes for studying black hole spacetimes and, indeed, proving without any dispute the Kerr identity of 
these massive objects. Estimates for the event rate of these inspirals are still uncertain, suggesting that LISA 

' could detect a number of events ranging from few, to hundreds or even over a thousand, during a 3 — 5 years mission. 
The detection of numerous events will potentially improve our understanding of massive black holes demography and 
population synthesis 0- For those rare cases where the black hole's environment is 'dirty' due to pronounced gas 
accretion (quasars, active galactic nuclei), an EMRI signal would provide valuable information on the structure of 
the disc itself (see Refs. Q^S)- to their nature, EMRI systems can be accurately studied using well-established 
black hole perturbation theory [Tol |. and after a substantial amount of work over the last decade or so, they can now 
be considered among the most well-understood sources of gravitational radiation. 

It is the purpose of this article (as part of a special Volume dedicated to gravitational radiation reaction) to review 
the main body of this work, in particular the computation of gravitational waveforms and the backreaction to the small 
body's geodesic motion. Deliberately, we only discuss results obtained by the 'flux-balance' method: the amounts of 
energy and angular momentum radiated away at infinity and the hole's horizon equal the amounts removed from the 
orbit. Unfortunately, this prescription only works for special classes of Kerr orbits: equatorial or circular-inclined. 
The orbits expected to occur in reality have no reason to belong to any of these cases. On the contrary, they are 
expected to be generic: with non-zero eccentricity and inclination. For these orbits, only two (energy and angular 
momentum along the black hole's spin axis) out of the three orbital 'constants' can be evolved using flux-balance. 
The third one, the so-called Carter constant, requires knowledge of the local gravitational self-force. The description 
of the inspiral in terms of this force involves a more sophisticated machinery than the simple flux-balance method, but 
at the end of the day, this higher precision will be required for constructing waveform templates for LISA. The self- 
force formalism is thoroughly reviewed in other contributions of the present Volume. In the meantime, flux-balance 
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calculations provide rigorous results for non-generic orbits, while certain approximations are still possible for generic 
orbits. In the future they will also serve as benchmarks for self-force results. We should warn the reader that due 
to space limitation and this Volume's particular subject, we do not discuss any data-analysis related issues. This 
omission should by no means reflect on the importance of data analysis for the LISA project; on the contrary, we 
believe that this crucial component deserves its own review exposition in the near future. 

The remaining of this article is organised as follows. In Sections Hll fc ITm we summarise previous well-known, as 
well as more recent material on Kerr geodesic motion. Next, in Section fllll we review the Teukolsky-Sasaki-Nakamura 
formalism on which all rigorous EMRI calculations are based. That Section includes presentation of frequency-domain 
results (equatorial-eccentric and circular-inclined orbits) as well as (more recent) time-domain calculations. Section llVl 
is devoted to the 'hybrid' method, an approximate scheme for generating generic inspirals and waveforms. The topic of 
Section |V] is somewhat unorthodox as it addresses the issue of possible non-Kcrr identity of the central massive body. 
A direct astrophysical application of EMRI calculations - estimates of 'recoil' velocities for merging massive black 
hole binaries - is presented in Sect ion IVll Our concluding discussion, with a brief summary and possible directions for 
future work can be found in Section IVIll A number of Appendices can be found at the end of the article. Throughout 
the article we use geometric units G = c = 1. 



II. GEODESIC MOTION IN THE KERR FIELD 



Not surprisingly, geodesic motion in Kerr spacetime is a well-studied subject since the discovery of the Kerr metric 
itself; an exhaustive discussion can be found in Chandrasekhar's classic textbook [lol |. Adopting a Boyer-Lindquist 
coordinate frame, the well known equations of motion are (with A denoting the affine parameter, related to proper 
time as X — t/ fi), 

E-^ = ±i?i/2^ R = T^-A\fi^r^ + (L,-aEf + Q], T = E (r^ + a^) - a 
dX 

= ±ei/2, Q = Q-cos''9[a'{ii^-E^)+Ll/sin^e] 



dX 



= - aE 



Lz \ aT 



dX V sin^dj ' A 

E-^ = -a{aEsm^9- L,) + {r^ +a^)T/A (1) 

where ^g^od = ] is the test-body's worldline and A — — 2Mr + a'^, Y. ~ + cos^ 6*. The 

energy E, angular momentum along the rotation axis L^, and the so-called Carter constant Q constitute the integrals 
of motion. They naturally emerge as separation constants in the solution of the Hamilton- Jacobi equation in Kerr 
spacetime (see Section III A|l . 

Generic Kerr geodesies can be parameterised by a triplet of constant orbital elements: the semi-latus rectum p, the 
eccentricity e, and the inclination angle t. For weak-field orbits {p 3> M) these elements reduce to familiar Keplerian 
notions. The elements p and e define the orbit's radial turning points, the apastron and periastron: 

= Tp = . (2) 

1 — e 1 + e 

The radial orbital period is defined as the time required for the body to move from apastron to periastron and 
back to periastron. 

In the strong field of a Kerr black hole, there are many ways that one could define an 'inclination angle' — for 
example, the turning points of the orbit's latitudinal motion, or the angle at which the small body crosses the equator 
as seen by distant observers. We use the following definition: 

Q = Ll tan^ t (3) 

This definition does not correspond to either of these examples, but is very convenient: it depends simply on orbital 
constants and has a useful intuitive description, suggesting that the Carter constant Q is essentially just the square of 
the angular momentum projected into the equatorial plane (this description is in fact exactly correct for Schwarzschild 
black holes; for non-zero spin it is not quite correct, but is good enough to be useful. We discuss this issue in more detail 
in Appendix IXji . The orbital elements can be written as functions of {E,Lz,Q}, and vice versa. These functions are 
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FIG. 1: Generic Kerr orbit with p — 2.5M, e = 0.5, t — 40° around a rapidly spinning hole a — 0.99M. The Figure displays a 
Boyer-Lindquist r — 9 coordinate slice of the orbit with the t/i-motion suppressed. 

obtained by solving for the radial motion's turning points, R{r) = 0. This equation is in general a quartic polynomial 
of r. Only the two largest roots are relevant for bound stable orbits, and they correspond to rp,ra- When l — one 
of the roots is trivially zero and we are left with a cubic equation. 

The actual turning points of the body's latitudinal motion are found by solving Q{9) = 0, which becomes, 

f3z^ - ( (3+ z + tan^ L = (4) 

\ COS"^ L J 

where z = cos^ 9 and j3 ^ a? (/i^ — £'^)/L^. From the two roots z_ < z+ the relevant one is the smallest and leads to 
the turning points On and 0s = tt — ^„ which are related to u as, 

en=i~^+0{a') (5) 

The latitudinal period Tg is defined as the time required for the body to travel from On to Os and back to On- 

Loosely speaking, Kerr bound orbits look similar to Keplerian ellipses especially when Tp > 15 — 20M. The main 
qualitative orbital features originate from the incommensurate nature of the three fundamental periods Tr,Tg,T^. In 
the weak-field Keplerian limit p » M these periods almost coincide leading to familiar elliptical orbits. At closer 
distances Tr — grows resulting in the famous periastron advance. In addition, provided that a 7^ 0, we also have 
Te 7^ which is responsible for the well known Lense-Thirring precession of the orbital plane. Both these effects 
become important in the strong- field region and any resemblance to a quasi-elliptical orbit is lost. An example of a 
strong-field generic orbit is given in Fig. ^ 

A point of the parameter space {p, e, l} does not necessarily correspond to a stable bound orbit. A separatrix surface 
Fs{p, e,L) = delimits stable from unstable (plunging) orbits. This surface is determined by requiring a double root 
for R{r) — 0, specifically r2 = Vp where r2 is the third largest root. For the case of Schwarzschild orbits the separatrix 
takes the simple form pT|. 

Fs=p- M{6 + 2e) = (6) 

The limits e = (circular orbits) and e = 1 (parabolic orbits) correspond to the familiar values = 6M and rp — AM 
for marginally stable and marginally bound orbits 10], 12]. For Kerr it is not generally possible to write down a simple 
expression as the functions E{p, e, l) etc. are quite cumbersome. For example, for equatorial orbits we obtain from 
r2=rp [a, 

F.^p" -{L,-aEf{l + e){i-e) = Q (7) 

but this is not as informative as IHJl. For the special case of prograde orbits around an extreme Kerr black hole, 
eqn. O simplifies io Fg = p — M{1 + e), which gives Vp = M for all eccentricities. For all other cases, the separatrix 
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FIG. 2: Left panel: Kerr equatorial separatrices, (from right to left): a — 0.5M (retrograde, dotted curve), a = (dashed 
curve), a — O.IM, 0.5M, 0.99M and a — M (dotted-dashed curve). For inclined orbits the behaviour of the separatrix is 
qualitatively the same. Right panel: number of revolutions per one radial period Tr for zoom-whirl orbits with e — 0.3 and 
black hole spin (from right to left), a = (dotted curve), a = O.IM, 0.5M, 0.99M and a = 0.999Af (dashed curve). Each of 
these curves terminate at a fixed distance 5p = O.OIA/ from the respective separatrix. 

is calculated numerically, see relevant plots in , 0| and |0| . In Fig. |2 we illustrate the location of the separatrix 
for equatorial Kerr orbits as a function of the black hole spin. When plotted on the p — e plane, the most prominent 
feature of the separatrix is that its location shifts towards smaller (larger) p with increasing black hole spin for 
prograde (retrograde) orbits. 

Of particular importance are orbits with p « Ps(e, i) (where ps solves Fs = ). These are called 'zoom-whirl' orbits, 
heuristically named after their particular properties p^ . A body in a zoom- whirl orbit wi ll sp end a considerable 
amount of its orbital 'life' close to the periastron. An approximation for Tr, a,s p —> ps, gives |llj.|l3|. 




which shows that the radial period will grow (and eventually diverge) as the separatrix is approached. In that region, 
the particle will trace a quasi-circular path before being reflected back to the apastron. Such behaviour will be 
particularly prominent for high eccentricity orbits: the particle will 'zoom in' from its apastron position, and perform 
a certain number of quasi-circular revolutions ('whirls') reaching the periastron (which should have a value close to 
Ps{e)/{1 + e)). Finally, the particle will be reflected and 'zoom out' towards the apastron again. Zoom-whirl orbits 
resemble a set of orbits known in the literature as homoclinic orbits [13. They occur in both Kerr and Schwarzschild 
geometries, and their potential significance for the detection of gravitational waves by space-based instruments was 
first pointed out some years ago by Cutler, Kennefick & Poisson ,1L] who concluded that the small number of whirls in 
the Schwarzschild case made the phenomenon less interesting for non-spinning black holes. The zoom-whirl behaviour 
is much more pronounced for prograde orbits around rapidly spinning holes, as illustrated in Fig. |21 

The number of revolutions plotted in Fig.|21is defined as = A0/27r, where A(j> = (j){t+Tr) — (l){t), the accumulated 
azimuthal periastron advance during one complete radial period. Using (jSJ, we find that AV exhibits the same 
logarithmic divergence as p — > Ps- According to Fig.|21 for small and moderate spins Nr stays close to the corresponding 
Schwarzschild value, but it grows rapidly as a ^ M, basically due to the intense 'frame-dragging' induced by the 
black hole's rotation in the very strong field region close to the horizon which can be reached by bodies in prograde 
orbits. Note that the spin dependence enters lO in the proportionality coefficient which we omitted. The overall 
behaviour can be understood as an extreme example of perihelion advance (as in the celebrated case of the planet 
Mercury) . 

In principle, as (jSJ suggests, the number of revolutions can be made arbitrarily large irrespective of the black hole 
spin, provided the body approaches sufficiently close to the separatrix. In practice however, this limit is immaterial 
as sufficiently close to the separatrix, radiation reaction makes a significant correction to the body's motion in each 
orbital period, and therefore we can no longer speak of slow, adiabatic, orbital evolution. The motion during the 
stage of separatrix crossing has been discussed in detail in Refs. 0, 0. 

In a realistic scenario, we should not expect to find (apart from chance cases where the particle enters a near- 
separatrix orbit as a result of its initial scattering) very high eccentricity zoom-whirl orbits, as it is well known that 
eccentricity tends to decrease under radiation reaction, for almost the entire inspiral. However, despite this decrease, 
a substantial amount of eccentricity will survive, in many cases, up to the point where the orbit is about to plunge 
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(the evolution of e during the inspiral is discussed in more detail later in this article). These orbits will probably 
become zoom-whirl orbits, especially when a rapidly spinning black hole is involved and the motion is prograde. Not 
surprisingly, the zoom-whirl behaviour will appear also when the orbit is non-equatorial (a good example is the orbit 
of Fig. Close to the separatrix we have S> le, and the body will spend a significant amount of time moving 
in a nearly circular-inclined fashion at r ~ rp. 



A. Hamilton-Jacobi theory in Kerr spacetime 



Despite our deep understanding of Kerr geodesic motion, certain aspects of it were studied only recently. For the 
most interesting case of generic orbits, an important issue concerns their periodicity properties. From this point of 
view, eqns. are somewhat misleading as they mix the r and 6 motions. As a consequence, it is not obvious how 
one could calculate orbital periods Tr,Te,T^ from these equations. However, as we are about to describe, the full 
separability of the Hamilton-Jacobi equation suggest that it is, after all, possible to decouple the t, r, 6, cj) motions by 
means of an action-angle canonical representation |l9j|.|l5||. In this way, the complete periodicity of Kerr motion is 
unveiled and rigorous fundamental orbital periods can be computed. 

The starting point is the 'super-Hamiltonian' for a test-body |23| , 



(9) 



where p°' = dx°'/dX is the four-momentum. Generating a canonical transformation {x°',pb} {/3"i76} with the 
requirement that the new Hamiltonian vanishes (and as a consequence the new coordinates and momenta are constants) 
leads to the Hamilton-Jacobi equation, 



(10) 



For a restricted class of coordinate systems (Boyer-Lindquist among them) this equation 'miraculously' admits full 
separation of variables (see j2ll for detailed analysis of this issue). The solution takes the form, 



S ^ -\- Et + L,(l) + Srir)+Sg{e) 



(11) 



with 



Sr{r) ^ I dr^ and Sg{e) 



de ve 



(12) 



Although the separation of the t and (f) coordinates is ensured by the stationary and axisymmetric nature of the Kerr 
spacetime, the separation of the remaining r,9 coordinates is not related to any obvious symmetry^. The Carter 
constant Q is related to the original r — 9 separation constant C a.s Q = — — a^E^. 
The equations of motion follow from 



dS 



p' 



'E- 
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t0 T 



(13) 



p0 ^ g4>4> L^-E g^* 



These are just the geodesic equations Similarly, the covariant momentum components are: pt 

±Vr/a, pg = ±Ve, p^ = L,. 



-E, Pr 



^ Unlike E and which are related to Killing vectors t", 0" associated with stationarity and axisymmetry (with E = —pa f^, Lz = Pa4'°'), 
the Carter constant is related with a rank-2 Killing tensor, Q = KabP"" |22| . 
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The tr ue p eriodic nature of geodesic motion in Kerr is revealed when we introduce a set of action-angle variables 
[m (see |l9j for background material). The actions are defined as (hereafter i = r,9, (f>) 



Jt^ f Pa dx" 



Jr = 2 dr Vr/a 

Je = 2 Jg" d9 VQ 



(14) 



For the full four-momentum vector we can choose, 



7a 



-E,Ji 



(15) 



which is to be considered as a function 7a = /a(— E, Lz, Q). A triplet of fundamental frequencies can be defined 
as, 



dU 
dJi 



(16) 



These are associated with motion with respect to r, 9, (j) coordinates. In (|16|l the Hamiltonian is assumed as a function 
of the canonical parameters {/3'^,7a}. In addition there is a fourth constant, associated with the t coordinate. 



m 

"dE 



(17) 



Unlike the familiar action-angle treatment of the Kepler problem jl9l | , it is not possible to integrate and subsequently 
invert eqns. H14|l to obtain 7i(7a). Nevertheless, as shown by Schmidt 15], it is still possible to obtain the desired 
derivatives dTi/d^a by employing a theorem on implicit functions for the function fa- If D[f] is the Jacobian matrix 
of / then the theorem states (provided that AetD[f] 7^ 0) that D[f] ■ D[f-^] = D[f] ■ {D[f])-'^ = I. Using the fact 
that H{—E,Ji) — — we can then obtain dTi./dE,dH/dJi as well as dQ/dE,dQ/dJi. Then, the fundamental 
frequencies (|lt)|) are found to be. 



K{k) 



ve = 



l(/.2_ ^2)1/2 ,^x{rp,ra) 



V4, 



K{k) Z{rp, rg) + Lz [n(zi, fc) - K{k)] Xjvp, r,) 
27rA 



2A ' 4A 
For the fourth frequency we find, 

dn K{k) Wirp, ra) + a^zlE [K{k) - E{k)] X{rp., r^) 



dE 



A 



(18) 



(19) 



The various integral quantities appearing in these formulae are listed in Appendix Note that the frequencies H18() 
are associated with the test body's proper time. The role of vt is to translate them to frequencies Sl^, fig, associated 
with Boyer-Lindquist time t io\. These are simply related as. 



The above frequencies are simplified by assuming small e and l. Then, 

1/2 



57^ — 



1 ± af^K 



1 =F 4arjK + 3^ 



P P 



21 1/2 
2 



(20) 

(21) 
(22) 
(23) 



where fix = yM/p^ is the Keplcrian frequency and the upper (lower) sign refers to prograde (retrograde) motion. 
Figure 121 shows how the corresponding periods Ti — 2n/Q.i vary with p for a near extreme hole a = 0.999M for both 
prograde and retrograde motion. The same qualitative behaviour persists for the general periods H2U|I . 
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FIG. 3: Fundamental Kerr periods Tr,Tg,T^ for small eccentricity, small inclination orbits. The black hole spin is a = 0.999M 
and both prograde (left panel) and retrograde orbits (right panel) are shown. Near each corresponding separatrix Tr oo 
(zoom- whirl behaviour) while Tg , retain a finite value. Only for the special case where rp ft; r+ (as in this Figure, left panel) , 
the ^-motion is 'locked' by the presence of the event horizon and Te — > oo too. For p ^ M all periods converge to the usual 
Keplerian period. The results shown here are typical even when the assumptions of small e, i are dropped jl^. 



The (constant) coordinates conjugate to 7a are given by, 



/3a = 



dS 



dW _ dW dak 
dja dak d-fa 



(24) 



where ak — (— E, L^, Q) and W{r, 6, (p, ak) — S — A/2. These equations define the Kerr angle variables Wa and 
provide the frequency interpretation for i^a- Expanding these we obtain 



dSr dSe 
dfi^ dijfl 



dSr dSg 

dE ^ dE 





F Mr A = 


— 2 Vr 


'dSr 


dSe' 
~ d^i^^ 


_^ 9Q 

dJr 


'dSr 

Jq^ 


dSe' 
dQ 








Pe- 


\-iygX ^ 




'dSr 

dy? 


dSg- 
dfi'^ 


dQ 
dJe 


'dSr 

Jq^ 


dSe' 
dQ 












-2v^ 


'dSr 


dSe' 


^dQ_ 
dJ^ 


'dSr 

_dQ ^ 


dSe 
dQ 


1 

+ 2^ 


dL^ 


dSe' 
dL,_ 



(25) 



By differentiating and combining these expressions one can arrive to the equations of motion For example, from 
the Wr^wg equations we get. 



= 1 



^2 R-y^f + a'cosHQ-^/'e^^^Kg 



= 1 



where the overdot stands for d/dX, and Ki = {dQ / dJi) / Vi. Combining these equations. 



{Kr - Kg) 



Q-^'H-R-^/'f 



(26) 
(27) 

(28) 



The function (T(r, 9) is subsequently specified by substitution in either H26|l or (|27() . We get a = {r"^ + a' cos' 6)^^ ^ 
In this way we have arrived at the expected r, 9 equations of motion. In a similar manner, the remaining t, t 
equations of motion follow from the Wt,w^ equations. 
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III. GRAVITATIONAL RADIATION FROM A TEST-BODY IN KERR SPACETIME 

A. Basic timescales 

Before venturing into our discussion of rigorous EMRI results it would be useful to provide some basic estimates 
for the timescales involved. The notion of adiabaticity in orbital evolution plays a central role in EMRI calculations. 
Quite simply this is a statement about the relative magnitudes of orbital periods and the radiation reaction timescale 
Trr. For the former we can use the Keplerian period, 



2t:M 



(1_ £2)3/2 \M 



(- 



3/2 



(29) 



which is close to the true periods Tr,Tg,T^ provided the orbit is not too close to the separatrix (see Fig. In this 
latter case, Tr becomes much longer than the other two and we should choose Torb = Tr. 

For the radiation reaction timescale we can choose among Tp = Te — e/\e\ or = The first two are 

comparable for p M while the third is always much longer (see Section llV)l . As the separatrix is approached Tp 



becomes the shortest, and it makes sense to choose Trr = Tp. From the weak-field expressions of Section Hvl we have, 



T 



RR 



T -JL 

' ~ \p\ 



M2 

— (1 -e2)3/2 



We can easily translate these timescales into human units. 



(30) 



T 



orb 



7.6 



(1 - e2)3/2 



RR 



M 
106 — 



Mepj 



(1 



2)3/2 



(31) 



where Afg — Af/lO^Af© and pe = p/6M. Clearly, Trr ^ Torb for the relevant (for LISA) mass ratios fi/M 
10"^ — 10~^. This condition only breaks down very close to the separatrix, p ps, as — > oo. 



B. Frequency domain calculations: The Teukolsky formalism 



For about three decades, the celebrated Teukolsky formalism has been the prime tool for the calculation of 
gravitational fluxes and waveforms from test-bodies moving in Schwarzschild and Kerr spacetimes. The eponymous 
equation describes the dynamics of linearised radiative perturbative fields in a Kerr background. In particular, instead 
of dealingdirectly with metric perturbations, the Teukolsky formalism considers perturbations on the Weyl curvature 
scalars [ly. The most relevant scalar V'4, is a result of the projection of the Weyl tensor on the null vectors n", 
fh^ which are members of a Newman- Penrose tetrad, that is, ^4 = —Cafj-ysn^rh^n'^rh^. The feature that makes this 
formalism so attractive for the EMRI problem is that the radiative fluxes (at inflnity and at the horizon) , as well as the 
two gravitational wave polarisations /ix can all be extracted from ■04. In a hypothetical absence of the Teukolsky 
equation, one would have to compute these quantities by solving the ten coupled metric perturbation equations. It 
is a great fortune that the Kerr spacetime, the unique spacetime describing astrophysical black holes, is a Petrov 
type-2? spacetime, which is the key property for arriving to decoupled perturbation equations for the Weyl scalars 
(see Ref. ^] for detailed analysis). As the Teukolsky formalism is a well-reviewed subject (see for example 25] for a 
detailed discussion in the context of EMRI), we shall only provide a summary and outline some basic results. 

The original Teukolsky 'master' perturbation equation is [2^, 



(r2 



-2s 



A 



2-2 

— - - a sm I 



M) 



A 



4 Afar 



90 



2s 



dtdcj) 
M{r^ - a2) 



a 
A 



1 



sm 



A-' — 
dr 



9-0 
'dt 



dr 



1 d 



sin 6* de 



sm 



de 



s\s coi^ t 



1] V = 47rEr 



(32) 



The relevant value for the spin-parameter is s = —2, which makes ip — p^^ip4, where p — (r — iacos6)^^ . The source 
term T is made out of projections of the given energy-momentum tensor, which for the present problem is the one 
for a point-particle. 



(33) 
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along the Newman-Penrose tetrad |23] (here u° = This is the part of the formaHsm where the small body's 

motion enters and is assumed to be known. Based on the adiabatic character of an EMRI system one is allowed to use 
the body's geodesic motion, unless the system is examined over time intervals comparable to the radiation reaction 
timescale. 

Eqn. (|32() is separable in the frequency domain by means of a decomposition, 
^4=/E/'^'^^^™(^'^)-2^"™(^)^"'"*''""^ and 4^ST-^y du;T,„(r,c^)_2 5F™We-*"*+^"'" (34) 

im im 

The function Rem satisfies the radial Teukolsky equation, 

^ T — ^ VrirjEim = Tim with Vt = r h SiLur + A, (35) 



dr \A dr 

where K = {r^ + a^) uj — ma and A = + a^w^ — 2amuj. The explicit form of Tgrn can be found in Ref. '2^. The 
spin-weighted angular functions -281^(9) (hereafter denoted simply as Sim) satisfy the following eigenvalue equation, 



1 


d 




sin^ 


)d9 





9 9 9 ^ iTT-^ . „ Am COS , 9„ „ ^ 

a''uj'^ cos^ ^ h Aau cos 9 H ^ 4 cot"^ 9-2 + Ei„ 

sin 9 sin 9 



Sim = 0. (36) 



These are normahsed as, 

d9sm9\Sim\'^ ^2tt (37) 



A particular solution of equation H35|l can be found in terms of two independent solutions of the homo- 

geneous equation, 

R^.n{r) = 1^ dr + I dr (38) 

where W the (constant) Wronskian Wl/S.^^/"^ R^^, A^^/^i?^]. Note that is vahd provided the source term Timir) 
does not extend to r cx3. Further discussion of this issue can be found in Section IlII Gl 

The solutions R^mt R^n ^"^^ chosen as to have, respectively, purely 'ingoing' behaviour at the horizon, and purely 
'outgoing' behaviour at infinity. Explicitly, 

rH ^ / A^e-^'^'''* R°° ^ I C'outgifcr. _^ ^2(Jin^-^kr, f^^, r^r+, ^on. 



^3 j^out^iuir, B^e" 



r e * lor r — )■ -f-oo 



where k — lu — ma/2Mr+, r+ = Af + (Af^ — the outer event horizon, and r, is the usual tortoise coordinate 

defined by dr^,/dr — {r"^ + a^) / A. From these expressions it follows that W = 2iwi3'". The solution H38|l describes 
ingoing waves at the horizon and outgoing waves at infinity as should be required on physical grounds. That is, 

A2 — ifcr» roo rp (^l\T>oo ( l\ 

^^'"(^ ^ " A^'" " ^wA^e-'^- (40) 

The source term T^m takes the form 0| , 

/oo 
e*"*-*™^(*) A2 Ai(r, 9) 5{r - r{t)) + {A2{r, 9) 6{r - r{t))} , -I- {A^ir, 9) 5{r - r{t))} „ (42) 
-oo ' ' J9=9(t) 

where Ai, A2, A3 are known functions of r, 9. The complex amplitudes defined in (I41II can then be written as, 

Zr^;" - ^-J^ / dtirm!^[rit),9it)] e-*-™^(*) (43) 
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where 



Rim A^-^—A, + ^-^A, 



(44) 

r=r(t),e=0{t) 



So far, all expressions listed in this Section are valid for any bound Kerr orbit. For the special case of equatorial 
orbits, we have 0{t) = -k/2 which makes functions of r{t) only. It is straightforward to show ^ that the 

quantities, 

= X°°'^(r(i)) e-'""!'^^*)-"**] (45) 

are periodic functions of time (with a period equal to the radial motion's period Tr). Consequently, they can be 
expanded in a Fourier series 

fc=-oo 

with ill- = 2TT/Tr. Plugging these in we arrive at 

^Zi" = ^Tnlk^i^ - '^rak) , (47) 

k—~oo 

where cumk = m -\- k^lr and 







case, I"^^ are functions of 9(t) only, and therefore periodic with period Tg. Then we arrive at. 



The last few steps can be repeated for another special case, that of circular-mclined orbits, i.e. r = p, e = 0. In this 

Emu) 

where now uJmk — fn ^0 + kQg. 

Inserting either H48|) or H49|) back in H34|) one can obtain the following expressions for ip4 at infinity and on the 
horizon, 

f Pir+)'^ J2imk ^Lk for r^r+ 
i'4^< _ (50) 
l'^ ^Etmki'Zk for +00, 

where 

^tt = ^^™(0)e-'-"-(*--)+""^ (51) 

Having at hand the Weyl scalar -04 we can immediately relate it to the two polarisation components hj^,hx of the 
transverse-traceless metric perturbation at r — + 00 , 

\ / ^„jj. mk 



The waveform's frequency spectrum is discrete, composed of the orbital harmonics cumk- 

Expressions for the gravitational wave energy and angular momentum fluxes at infinity can be derived using the 
Landau-Lifschitz pseudotensor [26| . 



dE^ 


00 


Itt ^ 


GW 


dL^ 


00 


~dt, 


GW 



1 /".,.„.„ o I / dh< \ ^ / dh-, ^ ^ 



IHt: W dt \ dt 



(53) 



— yd0d.sm^.^|^^ + — — ^ (54) 
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Time-averaging (over one orbital period or Tg) these expressions we get [27|. 

^GW - 2^ and Law - 2^ ^ 3 (55) 

emk ™fc imk 

The calculation of the respective fluxes at the black hole horizon is a subtle issue as expressions such as (|53|l . (|54f) 
are not available. Despite this difflculty, Teukolsky & Press were able to derive formulae for the horizon fluxes 
using the results of Hawking and Hartle j2^] on the change of the horizon surface area under a given gravitational 
perturbation. The end result is, 

Eg^^J^aerJ-^^ and ig^ = ^ (56) 

imk imk 

where a^mfc is a complex constant (its explicit form can be found in Ref. p^V 

For the case of generic orbits certain of the above steps cannot be taken as such. The reason is that the functions 
(|44|l are multiply-periodic and as a consequence their Fourier transform cannot be directly inverted. As recently 
discussed by Drasco & Hughes 29], this inversion is possible when the orbit is suitably re-parametrised [s^. A 
detailed analysis of generic EMRI will appear soon ■ 



C. The Sasaki-Nakamura formalism 



According to the previous discussion, the calculation of gravitational waveforms and fluxes boils down to the 
computation of the complex numbers and of the angular functions Sim- Before that, one needs to integrate the 

orbital equations and obtain {r{t), Bit), (/'(t)} and the frequencies fJ^, i^ei In addition, the homogeneous version of 

eqn. (|35|l has to be integrated, as i?^^ and their derivatives as well as the asymptotic amplitude B™ are required. 

Numerical integration of the Teukolsky equation itself proves to be somewhat problematic due to the long-range 
nature of the potential Vt- According to H39() . the B™ term decays at infinity much faster than the 5°"* term and 
can only be extracted with very low accuracy 0. In order to tackle this problem, Sasaki & Nakamura developed a 
formalism that maps the solutions of the Teukolsky equation to solutions of the so-called Sasaki-Nakamura equation 

^-F{r)^-U{r)X^O (57) 

The 'potentials' F{r), U{r) as well as the functions a{r), P{r), ri{r) appearing in subsequent formulae can be found in 
[25|. The solutions of this equation are related to the solutions of the Teukolsky equation via the differential operator. 



a 



A y (r2 + a2)i/2 Adr\ {r^ + a2)i/2 



(58) 



The key property of H57|) is that it encompasses a short-range potential. This can be demonstrated more easily if we 
shift to the function, 

Y{r) =r]^/'^{r)X{r) (59) 
Then, eqn. H57|l transforms into the Schrodinger-type equation. 



UY (60) 



with the short-range effective potential. 



4 2ry(r2 + a2)2 [ t] + J ) [k^ + 0{e'"^') for r* ^ -00, 

where c = (r+ — r_)/2M is a positive constant. As a consequence, eqn. (|57|l admits solutions, 

1 A"' e-'^''' + yl™' e*"*^* ' 1 Pe*'^'^* for r -f +00. ^ ' 
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FIG. 4: Examples of equatorial zoom-whirl waveforms (only h+ component shown, with the observation point located at 
(r, 6*) = (D,90°)) for the equatorial orbits p = 1.8M, e = 0.4, l = 0° (bottom panel) and p = 10.5Af, e = 0.5, i = 180° (top 
panel). The black hole spin is a = 0.99Af. 

The relation between the asymptotic amplitudes appearing in H39|) and (|62ll can be deduced from H58|l . 

S'" = --^A''^ (63) 

Unhke the Teukolsky equation, the Sasaki-Nakamura equation is ideal for numerical integration, as the 'ingo- 
ing/outgoing' components are of comparable magnitude at infinity and the horizon. We can then simply find the 
desired amplitude from H63|l . Similarly, knowledge of the wavefunction X{r) and its derivative at a given point 
immediately leads to the Teukolsky function Rtmir) and its derivative via the rule ()58|l . 



D. Equatorial orbits 

The Teukolsky-Sasaki-Nakamura formalism has been the standard machinery for computing gravitational waveforms 
and the accompanying back- reaction for test-bodies in equatorial orbits in Schwarzschild and Kerr spacetimes |llL IT^ 
IH) EE HE ^1 Oil ' Once the E'gw, ^gw fluxes are obtained at spatial infinity and the black hole horizon, 

then flux-balance dictates E = —Eq^/, — ^Lq-^. Since p — p{E, Lz), e = e(£', Lz) are known functions we find, 

p = iJjqi [L.^e E ~ E^.Lz] (64) 
e - H-^[-Lz,pE + E^pLz\ (65) 

where Foq = E^pL^^e - E^eLz,p- 

We first focus on the waveforms. For given values of p, e, provided rp is sufficiently large (something like > 20M), 
there is no significant difference between Schwarzschild and Kerr waveforms, since at such distances the black hole 
spin influences only marginally the orbital motion, and consequently the waveform. This similarity is no longer true 
in the strong field, especially when the hole is rapidly rotating. The waveform is significantly affected by the body's 
high velocity, intense frame-dragging and orbital plane precession, as well as by backscattered 'tail' radiation. 

Of particular importance are the 'zoom-whirl' orbits we discussed in Section ^] (more details can be found in 'l^). 
Provided that the eccentricity at this stage is not too small (which is the most likely scenario, according to the results 
of Section llVp . then the zoom- whirl behaviour should be taken into account for constructing waveform templates for 
LISA. A compact body in a zoom-whirl orbit will spend a considerable fraction of the orbital period in strong field 
regions and hence will radiate strongly. It turns out that a good fraction of the averaged flux is radiated during 
the motion near the periastron. The emitted waveform is accordingly shaped in a series of rapid 'quasi-circular' 
oscillations separated by relatively 'quiet' intervals. Once again, for values a ^ 0.9M the zoom-whirl effect will be 
enhanced (assuming a prograde orbit, see Fig.|21). Examples of zoom- whirl waveforms are shown in Fig. 01 generated 
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by the equatorial Kerr Teukolsky code of Ref. [I3. These waveforms correspond to two opposite extreme cases; 
prograde and retrograde orbits with p very close to Ps for a = 0.99M. In agreement with Fig. [3 the 'whirling' portion 
of the waveform is much more pronounced in the prograde case. Equally contrasting is the relative contribution of 
higher £— harmonics: the retrograde waveform is essentially composed hy £ < ^max ~ 4 — 5 while ^max ~ 12 — 14 for 
the prograde waveform. Another important aspect (not evident in Fig. ^ is the waveform's directionality. Viewed 
off the equatorial plane, the waveform suffers a substantial suppression in high frequency features. This is because 
the wave's higher multipole components (which are responsible for the small-scale structure) are mainly 'beamed' to 
directions close to the equatorial plane 13]. This effect is similar to classical synchrotron radiation and it can be easily 
'explained' noting that at the limit I = m '3> 1 the angular function (for a — 0) Sim ^ sin^ 6 and hence the radiated 
flux will exhibit asin2^0«e-^(«— an gular dependence near the equatorial plane (see Refs. [23| for early studies 
in gravitational synchrotron radiation). 

Next, we discuss the evolution of equatorial orbits under radiation reaction. This topic has been the subject of 
several studies ^ useful way to present the results is by plotting arrows (p, e) for each given orbit 

(p, e) of the p — e plane, as in Fig. [S] For all points, we have that p < which simply means that the orbit always 
shrinks, until the point where the separatrix is reached and orbit becomes a plunging one. On the other hand, the 
eccentricity exhibits a more interesting evolution. As suggested by weak-field computations [s^l (see eqn. H93I) of 
Section HV|) we typically find e < 0, that is, the orbit tends to circularise. This phenomenon of orbital circularization 
as a result of some form of dissipation is seen in many astrophysical situations, such as that of satellites whose orbits 
are decaying due to atmospheric friction. The reason is that the dissipating mechanism causes the body to 'drop' 
in its potential well, the usual geometry of which ensures that the orbital eccentricity decreases. However, for orbits 
in the vicinity of the separatrix we surprisingly find e > 0. This effect, which is confirmed by both numerical and 
analytical calculations, is intimately connected with the notion of a last stable orbit. As this orbit is approached at 
the end of the inspiral, the radial potential R (defined in eqns. ^ becomes shallower (as the minimum turns into a 
saddle point at plunge), and this tends to increase the eccentricity of the orbit. Shortly before plunge this mechanism 
overcomes the circularising tendency. 

Using Ht)5|) it is straightforward to show that e > sufficiently close to the separatrix. During this final stage of 
the inspiral the fluxes approximately obey the relation (following from eqns. (|55|l . (|56|) as a consequence of Tr ^ 
which makes iVmk ~ wfi^), 

E « fl^L, (66) 

which is characteristic of a circular orbit. This relation simply means that most of the energy and angular momentum 
is radiated during the body's 'whirling' at r « Tp, during which the radius hardly changes and there is a single 
dominant frequency and its harmonics. 
Using H66|l . we have near the separatrix 

e « H-qi [-E^, + L,^,] (67) 

We always have, < and consequently the sign of e is determined by the sign of the function [—E^e + ^(t> Lz,e] 
which encodes the body's geodesic motion close to the separatrix. For p ^ M this quantity is positive, leading to the 
expected e < result, but sufficiently close to the separatrix p ^ ps it takes negative values irrespective of e. Based 
on these observations, it follows that there must be a critical curve on the p — e plane where e = 0. We should expect 
to find this curve in the vicinity of the separatrix, as the e > effect is related to the notion of a last stable orbit. 
For example, for nearly circular orbits, this critical curve is located a.t Vr.pt = 6.68M for a = and as a — > Af it tends 
to 'coalesce' with the separatrix at Ps M, for prograde motion [l3l ISa. l36l|. In the same limit, their separation 
becomes maximal for retrograde motion, pcrit = 9.8M. As the numerical data of Fig. [S] demonstrate, this behaviour 
persists for arbitrary eccentricities: the eccentricity gain region shrinks (expands) for prograde (retrograde) motion, 
as the spin increases. 

The vectors in Fig. were generated using E = Eoo + E^i (and similarly for L^) with the option of setting 1^ = 1 
(total flux - dashed arrows) or ^ = (horizon flux neglected-solid arrows). No obvious difference between these two 
choices appear in the retrograde case, but there is a clear deviation in the prograde case, for near-separatrix orbits. 
Surprisingly, inclusion of the horizon fluxes slows down the inspiral. Responsible for this seemingly bizarre behaviour 
is the phenomenon of superradiance. We expand further on this issue in Section IlII Fl below. 



E. Non-equatorial orbits 

Gravitational radiation from a test-body in an circular-inclined orbit was first studied by Shibata 0|, more than 
ten years ago, by means of the Teukolsky-Sasaki-Nakamura formalism. Shibata not only computed fluxes E^ Lz and 



14 





FIG. 5: Radiative evolution of equatorial orbits, represented by arrows living on the p — e plane. The arrows have components 
(M/fi) [p, e]. The right (left) panel corresponds to prograde (retrograde) motion around a a = 0.99M Kerr hole. The separatrix 
and the e = critical curve are represented by the dashed and solid curve, respectively. Note the markedly different areas of 
the e > regions. 



waveforms but he also provided a first approximate expression for the Carter constant flux Q (more details are given 
in Section llV B|) . 

A few years later, Kennefick & Ori [45| were able to derive an exact expression for Q in terms of the components of 
gravitational self- force, see Appendix [CI for more details. The Kennefick-Ori formula relates the instantaneous fluxes 
{dE/dr, dL^/dr, dQ/dr}, 



-J- = Ti.E-j- +y-,L,-; 2 A^u,. 

dr dr dr 



The function H{r, E, L-^) is defined in Appendix lO 

For a circular orbit (r — const, Ur — 0) eqn. yields, 



(68) 



(69) 



for the averaged fluxes. As discussed by Kennefick & Ori, this is a key relation for the so-called ^circularity theorem'' 
|45 | .|25 | .|46 | : Kerr circular orbits remain circular under the influence of radiation reaction. From a practical point of 
view, eqn. H69|l is an exact expression which relates Q to the E, fluxes ~ quantities that can be calculated within 
the Teukolsky framework. Therefore, for circular orbits, the Q flux is within reach without having to resort to any 
self-force computation, or any kind of approximation. 

This important result was fully exploited by Hughes [T^ in his study of circular-inclined EMRIs. As in the previous 
case of equatorial orbits, strong fleld circular-inclined orbits generate the most interesting waveforms. A representative 
example (taken from Refs. 14], 47]) is shown in Fig. (top panel) for an orbit with parameters p = 7M, t = 62.43° 
around a rapidly spinning hole of a = 0.95M. The dominant feature in the waveform's pattern is the modulation 
induced by the body's 0— motion (Lense-Thirring precession). The same effect operates on a much longer timescale 
when we set a — 0.05M, keeping the same orbital parameters, see bottom panel of Fig.|H| This could be anticipated 
from the fact that the precession/modulation frequency is fig — il^ « zL2aM/r^, at leading order. 

Even more interesting are the results regarding the orbital evolution. In a fashion similar to Fig. |S1 the rates p, I 
are illustrated as vectors on the p — l plane, see Fig. The relevant expressions are, 



P = H-,^'A-L.^.E + E^,L,] 
^ = ^cTrc [-Lz,p E - E^p L^] 



(70) 
(71) 



with iJcirc = E^iLz^p — Lz^uE^p. It is always the case that p < (i.e. the orbits shrinks), but L shows a somewhat 
less monotonic behaviour. Provided that the black hole spin is a < 0.952Af we have Z > for all orbital parameters 
corresponding to stable circular orbits. Hence, the typical behaviour is a tendency of orbital angular momentum to 
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FIG. 6: Waveforms from Kerr circular-inclined orbits (only h+ component shown, with he observation point located at (r, 6) = 
(-D,90°)). The orbital parameters are p = 7M, l — 62.43° and the black hole spin is a = 0.95M (top panel) and a — O.OSAf 
(bottom panel). 



anti-align with the black hole spin, in agreement with the prediction of weak- field/slow- motion calculations j3S| |. As 
expected by symmetry considerations, both calculations predict L = when the black hole spin is turned-ofF. 

For rapidly spinning holes with a > 0.952M the above general rule is no longer true. Orbits that are 'horizon- 
skimming' [4^ (that is with p rj^) evolve under radiation reaction as to have i < |14| . This counter- intuitive 
behaviour can be explained by looking closer at the properties of circular orbits close to the horizon. Normally, 
we have the intuitive relation dLz/db < 0. However, horizon-skimming orbits violate this rule, i.e. they behave as 
dLz/di > 0. Since i — Lz (9Lz/9t)~^ and since Lz < always, it follows that i < for horizon-skimming orbits fl3|. 

The most striking feature in the radiative evolution of circular-inclined orbits is the small rate at which inclination 
changes, i.e. t/t <C p/p- In other words, during an inspiral, t appears to remain almost fixed (at least when e = 0) 
under radiation backreaction, especially when the body is not moving close to the horizon. The data presented in the 
following Section will make this statement more quantitative. 



F. Evolving the orbit 

We have so far discussed orbital evolution results in the form of {p, e) and {p, 1) vectors. Although quite instructive, 
they are, in a sense, 'frozen' in time. They do not convey a clear picture of an entire inspiral trajectory, which would 
start with some initial orbit and terminate at the point where the separatrix is crossed. Perhaps the simplest way 
to construct such inspiral trajectory is by pasting together a sequence of orbits. This should be a rather good first 
approximation, due to the adiabatic nature of the evolution. Such a program has been carried out by Hughes |47| 
for the case of circular-inclined orbits. A similar study for equatorial-eccentric orbits is in progress and should be 
completed in the near future. The orbital evolution algorithm is the following sequence of Eulerian steps: 

U+i = U + 5t (72) 
Pi+i = Pi +pSt (73) 
ij+i = Li + iSt (74) 

and the orbit's physical trajectory is described by — z^f,Q^[t, E(pi, Li), Lz{pi, Li), Q{pi, Li)]. The circularity theorem 
[2^,^^,^^ (Section lHI E|) ensures that e = throughout the inspiral. The selected time-step St should be consistent 
with the underlying assumption of adiabatic evolution, i.e. St <^ T^,Tg. 
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FIG. 7: Radiative evolution of Kerr circular-inclined orbits, represented by arrows on the p—L plane. Each arrow has components 
{M / ^)\p,i]. The panel on the left corresponds to a = 0.998M, while the panel on the right refers to a slower spinning hole 
a = 0.8M. The dotted curves represent the separatrix for each case. Note the i < arrows near the a — 0.998M separatrix, 
and the very long arrows in the a = 0.8M case, corresponding to orbits which are about to become dynamically unstable. 



In order to calculate the associated waveform one could in principle plug z'^(t) into the Teukolsky source term, 
and follow the steps discussed in Section [ill Bl The resulting waveform would be quasi-periodic, exhibiting a familiar 
'chirping' character during inspiral. 

In Ref. a somewhat different prescription is given for generating the inspiral waveform. The true periodic 
adiabatic waveform (as given by eqn. I52|l becomes 

h+-ih^ ^-Y 6-™^'-*+™^ (75) 

r ^ — ' uj , 

Imk "^f' 

where Z^mk, Sem{&) and ujmk = mfl^ + kflg have effectively become functions of time as they need to be 'updated' 
at each timestep ti . The above two prescriptions for generating inspiral Teukolsky waveforms are not fully equivalent 
to each other (in fact nobody has calculated a waveform using the former method) but, most likely, they would give 
similar results. 

Representative results from Hughes' study are shown in Figs.|Hl&|Hl As anticipated from the previous discussion, 
even for the strong field inspiral of Fig. |Slthe inclination angle remains almost fixed; the total accumulated change in 
L is of the order ~ 2°. The change in l is hardly noticeable for inspirals that terminate at 5M, as will be the case 
when the black hole is not rapidly rotating (see Fig. 5 in Ref. 47]). 

Perhaps, this is the most important result stemming from the study of circular inclined orbits. When combined with 
Ryan's weak-field result [i^ for i for the case of inclined and eccentric orbits (which also predicts that l/l <C p/p, e/e), 
it makes a strong case in favour of the assumption that l hardly changes even for generic Kerr inspirals. This idea 
was suggested by Cutler as a first stab at the problem of computing realistic generic inspirals, at least until self-force 
calculations are mature enough to be able to provide an accurate result for Q and i. In the meantime, the t = const 
rule is a practical approximation which could facilitate the computation of Teukolsky-based generic inspirals, along 
the lines of the calculation described here. The same rule is adopted in the so-called hybrid approximation of Kerr 
inspirals, which we expose in detail in Section llVI 

Another significant aspect of the inspirals of Fig. [S] is the importance of the fluxes at the horizon. According to the 
flux data from equatorial-eccentric 13] and circular-inclined orbits 14], the horizon fluxes can be as large as ~ 10% 
of the fluxes radiated at infinity. This extreme situation was illustrated in Fig. [S] for the evolution of equatorial orbits. 
Horizon fluxes of this magnitude can only be produced by bodies moving close to the horizon. In turn, such orbits are 
only possible around rapidly spinning holes, which also possess extended ergoregions. A fundamental property of an 
ergoregion is its ability to scatter back and amplify any incoming radiation field of frequency u < mflh = ma/2Mr^ 
(an excellent discussion on the many faces of superradiance and further references can be found in Ref. ^3)- Hence, 
it turns out that when the horizon fiuxes are significant they also represent superradiant back-scattered fields. In 
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FIG. 8: Teukolsky-based, Kerr circular-inclined inspirals for a a = 0.998M black hole for various initial inclination angles. The 
canonical values M — 10® Mq, ^ — IMq have been chosen. In generating these inspirals, the horizon fluxes were included (left 
panel) or ignored (right panel) . Each curve is labelled by the total inspiral time, starting at p = AM. 



effect, the orbit gains energy and angular momentum at the expense of the black hole's rotational energy'^. A perhaps 
more intuitive way to think of this effect is by considering the gravitational coupling of the small body with the tidal 
bulge that it induces on the horizon 113,1131 (described by the value of the Weyl scalar at the horizon). The 
bulge exerts a torque at the small body and tends to increase (decrease) its orbital frequency depending on the sign 
of $7^ — f2/i. This tidal coupling is not as exotic as it sounds, and is observed in other celestial systems such as our 
own planet and the Moon. 

The data provided in Fig.|Slfor the total duration of the inspiral (using canonical values /i = IMq^M — IQ^ Mq) 
suggest that the superradiance-induced inspiral delay can be ~ 5% of the total inspiral time. This effect is relevant 
only for rapidly spinning holes (a beautiful example of strong- field black hole physics that LISA will be able to probe), 
otherwise the horizon fluxes contribution quickly diminishes. 

Finally, in Fig. O we show the waveform corresponding to two different stages of the inspiral of Fig. |S1 (also taken 
from [43) with initial l = 40°. As the orbit becomes increasingly relativistic, the orbital plane precession frequency 
VLg — grows and strongly modulates the waveform. At the same time the monotonic increase of the rotational 
frequency U,^ shifts the waveform's frequency content to higher values. 



G. Parabolic and high eccentricity orbits 

The Teukolsky solution, eqn. H38I) . is not well-defined for orbits that extend to infinity, as the radial integrals over 
the source term diverge as r — > 00. In past calculations two methods have been used to remove this divergence. 
In the first method p39l |. the problem was dealt by carrying out integrations by parts, isolating the divergences in 
the resulting 'surface' terms. These terms were subsequently discarded. The second method involves the use of the 
Sasaki-Nakamura formalism, and in particular of the inhomogeneous equation, 

F{r) — -U{r)X = S{r) (76) 

The source term S{r) is related to the Teukolsky source Tim as [s^l, 

T^AWjr) ( . f, K\ d^W r^TUr) (. f, K\ 



^ In the early days of black hole perturbations, an interesting idea concerned the so-called 'floating-orbits' jgCj : bodies in strong-field Kerr 
orbits which would hover around the black hole without inspiralling, balancing their radiative losses at infinity by absorbing energy from 
the black hole. However it was soon realised that such conditions could never be realised. 
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FIG. 9: Teukolsky-based waveform (only h+ component shown, with the observation point located at {r,6) — (D,90°)), 
corresponding to the t ~ 40° inspiral of Fig. |H| 



Solving eqn. (|76|) via the standard Green function method and shifting back to the Teukolsky function we arrive at, 



Rimir -> oo) 



dr 



A 



(78) 



In contrast with Tim , the new source term decays as 



-5/2 



at infinity and therefore the above integral is finite. The 



price that one pays for this effective regularisation is the extra numerical integration of eqns. (j77(l . The above method 
has been used in a series of papers on equatorial parabohc/eccentric E3 ^^'^ circular-inclined orbits [s^ . 

The presence of divergent integrals in the Teukolsky solution (|38|l does not imply a failure of the Teukolsky equation 
itself. As explained by Poisson 53], writing the Teukolsky solution in the form (|38|l is legal as long as the integrals 
are convergent. If this is not true (as in the case of parabolic orbits), one needs to use the general expression. 



Rim{r) 



W 



A 



dr' 



Tim{r')Rfm 



ir') 



A2 



Rfmir) 



W 



B 



dr' 



T,m{r')Rrm{r') 



(79) 



where A, B,ri,r2 are constants. The relevant physical boundary conditions (purely outgoing/ingoing fields at in- 
finity/horizon respectively) can be imposed only after the integrals have been regularised. This can be achieved by 
integrations by parts |5.3| | . The resulting integrals are finite when we take the limits r2 — *■ c» and ri — > r^ and the 
divergent surface terms are then absorbed into the constants A, B which are eventually eliminated imposing the above 
boundary conditions. 

Eccentric orbits formally have well-defined integrals and therefore can be described using the original solution H38|l . 
However, for e close to unity the source inte gral s would grow (diverging at the e — > 1 limit). This remark is relevant 
for the Teukolsky equatorial codes of Refs. |llj.[l^ which are based on (|38|l . High-eccentricity tests have revealed 
that even for an e ~ 0.9 orbit there is a significant error induced by the growth in the source-integral. Therefore, 
for studying high-eccentricity and/or parabolic orbits one should either use the Sasaki-Nakamura solution (|78|l . or 
properly regularise the Teukolsky source integral ,3^ . 



H. Time-domain calculations 



Computing gravitational waveforms and fluxes via the solution of the Teukolsky equation in the frequency domain 
has been the 'traditional' approach to the problem for nearly thirty years, since the pioneering work of Davis et al [s^ 
and Detweiler '32]. This is understandable, as the full separability that the Teukolsky equation enjoys in the frequency 
domain essentially reduces the problem to the integration of two ordinary differential equations, eqns. (|35|l . (|36| ). 



19 



Working in the frequency domain indeed appears to be the appropriate computational tool for the EMRI problem 
since the frequency content of the radiation field is strictly harmonic (when geodesic motion is assumed). However, 
when the orbit is significantly different than circular equatorial, the emitted waveform is the sum of a large number of 
orbital frequency harmonics, especially for large eccentricity and/or inclination angle. In addition, strong-field motion 
requires the inclusion of several multipoles (£max ~ 10). Orbits with these properties are among those of relevance to 
LISA. For example, in the early stages of the inspiral all EMRI systems are expected to have e « 1 in which case any 
frequency domain code would require substantial computational resources. 

This has motivated some authors to approach the problem from a different perspective, that is, working directly in 
the time domain [s^ , [s^ , [s^] , [l^ ■ In the time domain, the Teukolsky equation is given by (|32|l . Only the Lp coordinate 
admits separation, hence we write 



-t-oo 



+ 00 



(80) 



Then H32I) can be rewritten in the form. 
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drl sin 6 89 
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- v{r, e) = t;,, 



(81) 



where $,„ = ^r,i\/{r'^ + a^) A^. The functions A, B, C,Tm can be found in Ref. |5q . 

Initial numerical investigations of the Teukolsky equation in the time domain was carried out by Krivan et al 
|59| for the homogeneous version of equation 181|l . Formulated as an initial value problem, these Teukolsky time- 
evolutions reproduced results for Kerr quasi-normal modes and late-time tails, already established by frequency 
domain calculations. 

Time evolutions of the inhomogeneous Teukolsky eqn. (|81|l were first performed by Lopez-Aleman, Khanna and 
PuUin H^ll^, for equatorial Kerr orbits. The code used in these studies is based on the original 2 + 1 homogeneous 
Teukolsky code of Ref. , appropriately modified to encompass a test-body source term. 

The time-domain representation of a test-body source term is a particularly delicate problem. The reason is the 
presence of (5— functions (and their derivatives) which are singular at the body's radial and angular instantaneous 
location (the (5((/) — (j){t)) function can be treated analytically). A first cut approximation is to smear out the S— 
functions by replacing them with narrow Gaussians (with width a of only few grid points) [55 | .|56 | . 
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(82) 



Typically, the radial grid is much more dense (~ few x 10^ points) than the 9— grid (~ 50 — 100 points ). 

In a time domain evolution, fluxes and waveforms are extracted by monitoring the field at some large distance (and 
close to the horizon for the horizon fluxes). The relevant flux formulae are [6l|. 
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(83) 
(84) 



As initial data the field is taken to be zero which is, strictly speaking, an inconsistent choice as it corresponds to a 
test-body appearing from nothing. This anomaly is manifested as an artificial, propagating, radiation burst at the 
beginning of the simulation (see Fig. 1101 below). Hence, one has to allow for sufficient time to pass in order for the 
computational grid to decontaminate . At the same time, the extent of the radial grid must be large enough as to 
prevent any spurious backscattering from reaching the 'observation' point, for a predetermined time window. 

The first results obtained in 55l._53| are quite promising as they generally agree within ~ 20% with the accurate 
frequency domain results of |l^,|43|, despite the crude representation (|82|) for the source term, and the two-dimensional 
nature of the problem. A sample of flux data, taken from Ref. (56(| , is provided in Table ^ In Fig. ^| we show a 
representative time domain waveform (also taken from Ref. H^) for the case of an equatorial eccentric orbit with 
p = 4.64M, e = 0.5. The black hole spin is a = 0.9M and the observation point is located at r = 150M, e = 7r/2. 

The second existing time-domain EMRI study is the one by Martel [Sg, and is focused on bound and parabolic 
Schwarzschild orbits. The separation of both (p and e dependence results in a 1 + 1 numerical scheme. This translates 
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FIG. 10: Time evolution of the real part of the Teukolsky function for m = 3 and = 0.01, for a test body in a 

p — 4.64M, e = 0.5 equatorial orbit, around a a = 0.9M Kerr black hole. The observation point is located at r = 150Af, 9 — ■k/2. 
Note the strong initial burst which is an artifact resulting from imposing initial data inconsistent with the presence of a test 
body. 



m mode 


E (Time domain) 


E (Frequency domain) 


Frac. diff. 


1 


2.9 X 10"" 


2.8 X 10"" 


0.036 


2 


1.4 X 10~"^ 


1.8 X 10""^ 


0.22 


3 


4.2 X 10"°® 


5.3 X 10"°** 


0.21 



TABLE I: Kerr time domain fiux data for the same equatorial orbit as in Figure nUl 



to a superior computational performance and smaller discretization-induced errors, as compared to the 2 + 1 Kerr 
calculation. 

Instead of the Teukolsky equation, Martel is evolving the inhomogeneous Zerilli and Regge- Wheeler equations which 
describe polar (e ven) and axial (odd) metric perturbation components, respectively. Both these equations have the 
well-known form 



yen) 



V{r)'lp£m = Stm (85) 



df^ dri 
with a source term, 

S,^ = Gi{t, r) 6{r - r{t)) + G2{t, r) 5' {r - r(<)) (86) 

The functions Gi,G2 can be found in |^. The numerical method used in evolving eqns. (p^ is described in jH3|, |H3 
and includes a state-of-art discrete representation of the source term's (5— function and its derivative. Martel's results 
are in excellent agreement (at the level of < 1%, see Tabled) with the corresponding frequency domain results, and 
not surprisingly, significantly more accurate compared to the Kerr time domain results. 

Apart from several tables of flux data, Ref. [s^ also presents several time-domain waveforms for bound and parabolic 
orbits (e < 1); these share the same qualitative features with the waveform of Fig. EH hence we do not present them 
here. 

Computing averaged fluxes in the time domain involves certain error, unlike in the frequency domain where the 
averaging is performed in an exact analytic manner. For circular and small eccentricity orbits the orbital periods are 
short enough (in comparison with the total duration of the time evolution) to allow for averaging over several orbits. 
In such cases the induced error is small. On the other hand when the eccentricity is substantial the radial orbital 
period becomes considerably longer, resulting in a growth of the averaging error. This situation is evident in the 
sample data of Table ITU Hence, high eccentricity orbits require long time evolutions in order to obtain high precision 
flux data. Meanwhile, as pointed out in "SS*], the waveform calculation has no such problems, as no time averaging is 
involved. 
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orbit 


flux 


Ireq. domain 


time domain 


frac. diff. 


V 


= 7.50478 


is 


3.1680 X 10""' 


3.1770 X 10"* 


0.3% 


e 


= 0.188917 




5.9656 X 10"^ 


5.9329 X 10"^ 


0.5% 


V 


= 8.75455 


is 


2.1008 X 10"'' 


2.1484 X 10"'' 


2.3% 


e 


= 0.764124 




2.7503 X 10"^ 


2.7932 X 10"^ 


1.6% 



TABLE II: Comparison of averaged fluxes for eccentric Schwarzschild orbits obtained in the time domain [H^. and the frequency 
domain [Tll |. 



As we mentioned, the main goal of the time domain approach is to be able to deal (primarily) with highly eccentric 
orbits, for which the standard frequency domain formalism becomes computationally expensive. For the Schwarzschild 
case there is no doubt that, indeed, the time domain approach is superior, partially because time evolving 1 + 1 wave 
equations is easily handled by modern computers and partially because there is an accurate representation for the 
singular test-body source term |63 | . For the more relevant Kerr problem, frequency domain calculations still lead the 
race, at least until Kerr time domain calculations reach the level of accuracy of the corresponding Schwarzchild codes. 
Indeed, there are ongoing programs that aim towards that direction (for example ,63.] ) . interfacing at the same time 
with self- force computations. 



IV. THE HYBRID APPROXIMATION 



In Section HTlFl we discribed a simple algorithm for constructing the orbital inspiral within the framework of the 
Teukolsky formalism. It is clear that this kind of calculation is computationally expensive as one needs to compute 
{p, e, Z} at each timestep. Using template banks of such Teukolsky-based inspiral waveforms seems as an unnecessary 
complication for current approximate LISA data analysis studies. Instead, it would be highly desirable to have at hand 
an alternative (albeit approximate) 'quick and dirty' scheme for computing the orbital evolution and the associated 
gravitational signal. 

Such a scheme is the 'hybrid approximation' js^l which is the topic of this Section. Essentially, it is a resurrection 
and extension of the older 'semi-relativistic' approximation by Rufhni and Sasaki |65j| . and is based on the following 
two-level idea: 

At the first level, in calculating the radiative inspiral we replace the fluxes {i?, i^} with approximate Post-Newtonian 
analytical expressions {E^^ , L^^} (hence removing the most computationally demanding component of the Teukol- 
sky calculation). For inclined orbits one has to handle the extra flux Q. This is not a problem for the special case 
of circular-inclined orbits where, as we have discussed, Q is fully of the other two fluxes. However, for generic orbits, 
the calculation of Q is elusive without resorting to the gravitational self-force. It is still possible, however, to come 
up with certain approximations for Q. We further expand on this issue below. 

In calculating the rates {p, e} we use the exact functional dependence on the fluxes, i.e. 

p = H-\bpE + CpL,+dpQ) ^ p = H-\bpEP^ +CpL^^ +dpQ) 

e = H-\beE + c^L^+deQ) => e ^ H-\be E^^ + L^^ + Q) (87) 

Note that the quantities H,bp_e, Cp^e, dp^^ encode the exact geodesic motion and are just combinations of derivatives 
of {E{p, e, i), Lz{p, e, l), Q^p, e, /,)}, see Appendix IdI For i we can write the simpler expression. 
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The inspiral {p(t), e(t), i(t)} in the space of orbital elements is then constructed by an Eulerian algorithm, similar 
to the one given by eqns. 1741) . The inspiral in the physical coordinate space {r, 0, 0} results from the combination 
of the full general relativistic geodesic equations (which are numerically integrated) and {p(t), e(t), i(t)}. With this 
prescription, the small body's world-line is described by the four-vector, 

^'^W = ^geod[i;^^W,i.W,QW] (89) 

where E{t) = E\p{t), e{t), t(t)], etc. 

At the second level of the hybrid scheme, we can compute 'kludge' waveforms by assuming that the small body is 
moving along the trajectory H89|l (with the option of taking into account radiative backreaction or not) like a 'bead 
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on a wire' in a fictitious flat spacetime, pretending that Boyer-Lindquist coordinates are true spherical coordinates. 
Effectively, calculating the waveform becomes a problem of wave dynamics for metric perturbations in flat-spacetime 
(see [23 for in depth discussion). Strictly speaking, this prescription is not part of a self-consistent approximation 
framework; it is rather more like a 'black box' designed to produce waveforms with accurate phasing. This is the most 
important property of a template waveform, and is crucially sensitive to the precision of the orbital motion. Indeed, 
the resulting waveforms agree well with the accurate Teukolsky waveforms. More details on the actual procedure and 
examples of kludge waveforms are given in the following Section. 



A. Equatorial inspirals and waveforms 



We first apply the hybrid method to equatorial orbits. Since Q — Q — 1 — l = the rates (H^Tjl reduce to eqns. Ipi^ . 
A good first choice for the weak-field fluxes are the expressions derived by Ryan for generic Kerr orbits. These 
are leading-order Post-Newtonian (PN), amended with the leading order spin term (which appears at 1.5 PN order), 
and correspond to a 'Keplerian plus leading spin effect' description of the orbital motion. In fact, at this level of 
accuracy, Ryan was able to use the leading-order self force (discussed in and derive a leading order result for Q 
(see Section riVB|) . Rewriting Ryan's fluxes in terms of our orbital elements yields (here we keep l arbitrary as the 
same fluxes will be used for non-equatorial inspirals). 
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(90) 
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The functions /i(e), /2(e), 51(e), 52(e), 03(e) are given in Appendix IeI In the limit a = 0, eqns. (|90|l and (|91ll reduce 
to the classic Peters-Mathews formulae |37l |. 

It is instructive to compare the hybrid inspirals with the inspirals generated by approximating the entire p, e, L 
formulae. At leading order we have. 
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These can be combined to give, 



P(e) = Pi 
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(92) 
(93) 

(94) 



where pi and Ci are arbitrary initial values. 

A sample of representative results for astrophysically relevant initial parameters is shown in Fig. 1111 We compare the 
leading-order trajectories (|94|l with the trajectories predicted by the hybrid scheme. The inspirals we show correspond 
to both prograde and retrograde orbits, for black hole spins a = 0.5M and a = 0.9M. In the weak-field region p'^ M 
the hybrid and the leading-order calculations agree well, as expected. Crucial differences between the two methods 
become apparent in the strong-field regime, during the inspiral's final stages (this is also the epoch most relevant for 
LISA'S observations). As is clear from eqn. the leading-order inspiral trajectory exhibits constantly decre asing 
eccentricity. This is not what the rigorous strong-field calculations, both numerical and analytical ITsL l35l l36j| . 
predict. 



As we discussed in Section Fill Dl there exists a region near the separatrix of stable/unstable orbits where e 
reverses sign: the eccentricity should grow near the separatrix. This feature is indeed present in the hybrid inspirals. 
Moreover, the location of the critical points where e = is in good agreement (at the order of few percent) with 
the numerical results of Refs. 1^,0]; see Ref. for actual data. Detailed comparison of the {p, e} values against 
accurate Teukolsky results establishes the superiority of the hybrid scheme over the leading order formulae (|94|) . The 
same conclusion holds even if we consider higher order PN expansions for p, e. A sample of representative data, taken 
from Ref. 64] is provided in Table ITTll 

Naturally, the hybrid approximation has its limitations. This can been seen, for example, in Fig. 1111 Three of 
the four cases shown appear 'good' in the sense that the trajectories appear to agree reasonably well with what we 
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a/M 


p/M e 


Calculation 


{M/ii)p 


Frac. diff. in p 


(M2/^)e 


Frac. diff. in e 





7.505 0.189 


Numerical 
Hybrid 
Leading order 


-7.475 X 10^2 
-6.859 X 10"^ 
-2.957 X 10"^ 


0.0824 
0.6044 


-1.967 X 10-3 
-1.291 X 10-3 
-1.159 X 10-3 


0.3434 
0.4108 





6.9 0.4 


Numerical 
Hybrid 
Leading order 


-4.240 X 10"^ 
-3.056 X 10"^ 
-3.420 X 10"^ 


0.2792 
0.9193 


+1.047 X 10-2 
+1.506 X 10-2 
-2.929 X 10-3 


-0.4384 
1.2797 


0.5 


6.5 0.4 


Numerical 
Hybrid 
Leading order 


-5.999 X 10-2 
-4.606 X 10^2 
-4.091 X 10-2 


0.2322 
0.3181 


-5.155 X 10-3 
-3.356 X 10-3 
-3.719 X 10-3 


0.3490 
0.2786 


0.5 


15 0.4 


Numerical 
Hybrid 
Leading order 


-3.371 X 10"=^ 
-3.358 X 10-3 
-3.328 X 10-3 


0.0039 
0.0128 


-1.341 X IQ-"* 
-1.334 X 10-"* 
-1.311 X 10-"* 


0.0052 
0.0224 


0.5 


4.8 0.3 


Numerical 
Hybrid 
Leading order 


-6.354 X 10-^ 
-4.858 X 10-^ 
-4.849 X 10-2 


0.2354 
0.9237 


+1.369 X 10-2 
+3.519 X 10-2 
-4.432 X 10-3 


-1.5705 
1.3237 


0.9 


5 0.4 


Numerical 
Hybrid 
Leading order 


-7.507 X 10-2 
-4.617 X 10-2 
-2.698 X 10-3 


0.3850 
0.9641 


-9.266 X 10-3 
-1.965 X 10-3 
-9.732 X 10-^ 


0.7879 
0.9895 


-0.99 


10.5 0.4 


Numerical 
Hybrid 
Leading order 


-7.516 X 10-2 
-5.506 X 10-2 
-9.704 X 10-3 


0.2674 
0.8709 


-5.312 X 10-* 
-5.223 X 10-"* 
-5.461 X 10-"* 


0.0168 
0.0281 



TABLE HL Comparing the rates p, e for some Kerr equatorial-eccentric orbits. The fifth column in this table shows p; the 
seventh column shows e. Within each section of the table, the first row of columns five and seven contains accurate numerical 
data from [TH . H^ . the second row shows data using the hybrid scheme and the third row shows data using leading-order 
results. The sixth and eighth columns show the fractional differences between the two approximation schemes and the accurate 
numerical results. In all cases but one, the hybrid approximation is closer to the accurate numerical calculation, sometimes 
substantially so. 



expect based on strong- field numerical analyses |l3j |. However, things clearly go wrong in the fourth case (a = 0.9Af, 
prograde); both the eccentricity growth near the separatrix and the distance of the critical curve e = from the 
separatrix are excessive. This behaviour is not surprising as prograde orbits of rapidly rotating black holes reach 
rather deep into the black hole's strong field where the weak-field fluxes (|9()(l and cannot be trusted. As a rule 
of thumb, and based on the data of Ref. l63| the hybrid method fails when Vp < 5M. This effectively constrains the 
black hole spin to a 0.5A/ for prograde motion. For retrograde orbits the hybrid approximation is favoured, since 
Tp never comes close to the horizon, regardless of the spin. 

An important prediction of the hybrid approximation is that for equatorial orbits the residual eccentricity prior to 
plunge should be substantial, in strong contrast to the prediction of the leading order formula (|94|) . In many cases, the 
leading order results predict that the orbit will actually circularise prior to plunge. Because the frequency structure 
of a nearly circular inspiral is rather different from that of an inspiral with substantial eccentricity, these results have 
strong implications for the waveform models to be used in LISA's data analysis. 

It is possible to further improve the performance of the hybrid scheme by simply using more accurate fluxes for 
E and Lz- As pointed out, Ryan's fluxes include only the leading (1.5PN) spin-dependent term; at the same time 
they are fully accurate with respect to eccentricity. Higher order (up to 2.5PN order) fluxes have been derived by 
Tagoshi [b^ , under the assumption of small eccentricity (see Ref. j23 for a complete collection of PN fluxes in the 
extreme mass ratio limit). In practice this is not a serious restriction as for the majority of inspirals the eccentricity 
is ^ 0.5 at the stage where the orbit resides in the black hole's strong field regime. At the same time, during the 
early stages of the inspiral we have rp ^ M and the higher PN corrections are small, irrespective of the eccentricity. 
However, the e— dependent factor in the leading PN part of the fluxes is important, but fortunately, this factor is 
fully known from the work of Peter & Mathews (33 • 

Hence, the optimal choice for the hybrid approximation would seem to be a combination of Ryan's and Tagoshi's 
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FIG. 11: Examples of hybrid Kerr equatorial inspirals, into a hole with spin a — 0.5Af (left panel) and into a hole with a — 0.9M 
(right panel). In each panel, the top half shows prograde inspirals and the bottom retrograde inspirals. In each set, the dashed 
line represents the separatrix separating stable from unstable orbits. The hybrid approximation is used to radiatively evolve 
orbits with initial parameters {pi,ei) = (20Af, 0.4), (20M, 0.99), (14M, 0.99), and (lOM, 0.99). The hybrid inspiral trajectories 
are shown as the solid lines in each plot (The final set is not included in the retrograde inspirals since the initial conditions are 
not stable in those cases.) The dotted trajectories in each plot shows the leading-order predictions generated using eqn. 1941 1. 
Note the significant qualitative difference between the two calculations at the vicinity of each separatrix. Note also the extremely 
large growth in eccentricity seen in the prograde inspirals for a = 0.9M just before reaching the separatrix. Comparison with 
accurate strong-field numerical results shows that this growth is excessive. Note that the time dependence of the inspiral is 
suppressed in this figure: most time is actually spent at large p. The total duration of an inspiral scales with M'^ / fj,. The 
shape of a curve, however, does not depend on this ratio, provided that the mass ratio is extreme; these curves are universal 
for < M. 



fluxes. In terms of our orbital elements, this union gives (full details will appear in |67l|). 




where the various additional coefficients are listed in Appendix IeI These fluxes include terms up to 2.5 PN order. 
In fact, some testing reveals that it is the 2PN fluxes that actually have the superior performance. This well-known 
result is related to the asymptotic nature of the PN expansion, see 25]. 

Another crucial improvement concerns the behaviour of the inspirals in the e ~ region. Although not particularly 
visible in the examples of Fig. 1111 a closer examination of that portion of the parameter space reveals that the hybrid 
approximation predicts a growth of e instead of the correct e ~ e behaviour _25l. ^35l._36j.,^45.] for nearly circular 
orbits. It is possible to rectify this problem by adding certain terms in the fluxes |2|,[63- Since the same pathology 
also appears when i ^ 0, we provide more details about this fix in our discussion on generic inspirals. Modifying 
accordingly the fluxes (|95|l . ()96f) . we can generate improved hybrid equatorial inspirals, as the ones shown in Figure [T^ 
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a = . 9 M prograde 




FIG. 12: Improved equatorial inspirals generated with the help of the new fluxes l|95|l . I|96|l additionally modified according to 
the low eccentricity fix discussed in the text. These inspirals are to be compared with the corresponding a = 0.9M prograde 
inspirals of Fig. 1111 The dashed curve represents the separatrix. 



We have chosen to display a = 0.9M prograde inspirals, a case where the original hybrid calculation was clearly 
inaccurate. 

The new inspirals look clearly more 'sane' all the way to the separatrix, devoid of the excessive eccentricity growth 
seen in Fig. 1111 Note that the residual eccentricity prior to plunge is e 0.1 not much different from what the leading 
order inspirals would predict. As we have mentioned elsewhere, for prograde orbits around rapidly spinning black 
holes, the e > region is small and the e = curve is located near the separatrix, hence eccentricity decreases almost 
for the entire inspiral. 

Having discussed the first level of the hybrid approximation (construction of inspirals) we move on to the second 
level, the generation of kludge waveforms 68]. We only present examples where backreaction is ignored, in which case 
it is straightforward to compare with the available Teukolsky-based equatorial waveforms of Ref. 13]. The starting 
point is the quadrupole formula for the flat-spacetime, trace-free metric perturbation hat — hat — ^riabh, 

y^ = -P''{t^r) where P'=(i) = / dF'x'^x"^ r°°(t, x'), (i, j. A: = {1, 2, 3}) 



(97) 

is the binary's quadrupole moment and T^'' is the small body's Minkowskian stress-energy tensor. All integrations 
are performed over the source. According to the hybrid prescription, 

z°-{t) = [t,r{t),e{t),(l){t)]Bh =^ z'^it) = [t,x^ = r{t) sme{t) cos (j){t) , x'^ r{t) sm0{t) sin ^(t) , = r(i) cos 6'(i)] (98) 

From h^'^ one can extract the transverse-traceless components /ix in a standard way po| . 

It is straightforward to enhance the accuracy of the hybrid waveforms by including additional source multipole 
moments. For example we can use the following quadrupolc-octupolc formula |27j].|69|. 



(99) 



where 



3"^''= dV'x''x"'T°''{t,yL') and AP^'' ^ dV'x'^x'^ x"'T°°{t,x') 



(100) 



are the current-quadrupole and mass-octupolc moments and = x /r. Indeed, it is possible to take into account all 
multipole moments by employing the 'fast' solution found by Press [70|, but it turns out that the simpler quadrupole- 
octupole waveforms show (in most cases) little deviation from the Press waveforms. In Fig. 1131 we present a small 
sample of kludge waveforms (for Kerr equatorial orbits) and compare them with the corresponding Teukolsky-based 
waveforms. Typically, the agreement between the two calculations is surprisingly good, despite the simple-minded 
construction of the hybrid waveforms. The comparison can be quantified in terms of the overlap function between the 
waveforms; in the present case the resulting overlaps are very close to unity (see Ref. for more details). The use 
of exact geodesic motion ensures high precision phasing for the waves, even in strong-field conditions. The agreement 
in amplitude is less impressive (mainly due to the negligence of higher multipole contributions and backscattering 
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FIG. 13: Sample of hybrid 'kludge' waveforms for the same Kerr equatorial orbits as in Fig.0]for a a = 0.99M Kerr hole. For 
the retrograde orbit p = 10. 5M, e = 0.5 (top panel) there is good agreement in both phase and amplitude with the analogue 
Teukolsky waveform. For the strong-field prograde orbit (bottom panel) p = 1.8M, e = 0.4 the two waveforms agree well only 
with respect to phase. Higher multipoles and backscattering, two factors not taken into account in the hybrid calculation, play 
an important role in this example. 

effects) but this has no significant impact on the data-analysis performance of the hybrid waveforms. We need to point 
out that the high degree of phase accuracy will not be preserved when the radiative trajectory H89|l is used (assuming 
the waveform is monitored for sufficiently long timescale ~ Trr), nevertheless, hybrid waveforms still remain valuable 
tools. A detailed quantitative discussion on the accuracy /reliability of hybrid waveforms will appear in Ref. |6l| . 

The punchline of this discussion is that hybrid waveforms have two key properties: they (i) can be very easily 
generated, for any Kerr orbit, and (ii) compare well to the existing rigorous Teukolsky-based waveforms. Hence, 
without much doubt, they can be used to produce a 'leading-order' waveform template bank for LISA. Indeed, the 
most recent estimates on the expected EMRI event rate were computed with the aid of hybrid waveforms . 



B. The elusive Q flux 



Attempting to apply the hybrid approximation to inclined orbits, one immediately faces the issue of calculating the 
Carter constant flux Q. Unfortunately, unlike E and this quantity cannot be inferred by monitoring the radiation 
field at infinity. Instead, in order to compute Q, one needs knowledge of the local self-force (see eqn. (|108|) below). 
In fact, one of the main motivations behind the self-force program has been the calculation of this elusive flux. As 
we already mentioned, Ryan |46j was able to derive the following leading-order result for Q making use of the known 
leading-order self force 



Qr 



64 „ [mV , 



2^3/2 sin^ 



51(e) 



a 

M 



M 



3/2 



/9(e) COSi 



(101) 



with /9(e) = 85/8 + (211/8)e2 -|- (517/64)6-^. Currently, this is the only available PN result for Q. A different 
approximation emerged from the study of radiation reaction for circular- inclined Kerr orbits [l3 | , which we reviewed 
in Section IIII El For this special orbital family, Q can be expressed as a linear combination of E, in an exact 
sense; no self-force is required. One of the main results of Ref. is that the inclination angle t remains almost 
fixed during the radiative inspiral, even for strong-field motion. In Appendix (taken from Ref. 64] ) we provide an 
intuitive explanation as to why the t = const rule should be a good approximation for all generic Kerr orbits. 
It follows from (|88|l that this rule yields. 



_ 20 . 

sph — " Y 



(102) 
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It is important to note that this 'spherical' rule is in fact exact in a spherically symmetric spacetime, like Schwarzschild. 
In such spacetime, the Carter constant is nothing more than the projection of the angular momentum on the equatorial 
plane, i.e. Qsph — L^ + Ly. Then it can be shown that an inspiral in this spacetime would proceed at exactly constant 
inclination angle, gravitational waves remove the amounts of L^; and Ly required to hold t fixed. 
Expanding now H102|l at the same order as Ryan's expression H101|) we get. 



5 VP/ cos i 



32 (e) cosi + — ( — ) {52(e) - /9(e) cos l} 



(103) 



At leading order this coincides with eqn. (|l()l|l . but there is disagreement in the spin-dependent term. Since by 
construction H101|) is fully accurate at 0{a) order, we conclude that the corresponding order term in (|1U3|I is incomplete. 
The missing piece is, 

Sq ^ Qn - Qsph = ^A^^ ^ ^ff2(e)(l - e^/' (104) 

5 \ P J M cosi 

This piece is also the leading 'aspherical' contribution to Q. As such, we would expect, 

A A ,2(5-64 ^ f mV''^ a sin^t 
Lz 5 \P J A-i cosi 

to be an improvement over the spherical formula 1)102(1 . see Ref. Instead of keeping l constant, this new formula 
predicts an increase, 

L,Q~2QL, 32 fi^ ( a 



2VQ{Q + Ll) 5 {Q + LD 



(^)'^sin.(l-e^)3/2g2(e)>0 (106) 



This is a welcome feature, in qualitative agreement with the accurate Teukolsky results of Section Fill El 

Before moving on, we should mention another weak-field approximation for Q, the first in chronological order, which 
was proposed by Shibata 01 for the case of Kerr circular-inclined orbits (before the formulation of the circularity 
theorem). In terms of L± = L^, ± iLy and their averaged fluxes L± , 

Qk{L+)L+ + {LJj L_ - 2a^ E{cos 9) E (107) 

where the time average is taken over a period Tq. This expression (which nowdays is obsolete), together with exact 
Teukolsky-based fluxes E,Lz allowed Shibata to compute the backreaction to the orbit (in addition to waveforms), 
with results similar to Hughes' ^14]. 

At this point we have at hand three different approximate expressions for Q, cqns. H101|) . H102|) . H105|) which can be 
incorporated in the hybrid scheme and construct generic Kerr inspirals fSection IIV C|) . 

Further insight into Q is provided by considering the exact formula for the instantaneous evolution of the Carter 
constant as given by Kenncflck & Ori [43, see Appendix [HI This time, Q is expressed in the form (|C2|) which leads 
to, 

dO „ dE „ dLz „ 9 „ 9 „ dE „ dLz ^ 

^ +g ^ +2pgFe ^ -2a^Ecos^9 — + 2LzCot^9 —^ + 2peFg (108) 

dr rfr ' ~ dr dr dr 

This is of course equivalent to expression H68|) . It is an interesting exercise to derive from eqn. H108|) the spherical 
rule (|102|l . The following argument is taken from Ref. |63|. In the a — > limit, the first term vanishes, the Carter 
constant reduces to the (square of the) projection of angular momentum on the equatorial plane, Q = — and 
the inclination angle l becomes the true inclination of the orbital plane, since tan^ t = Q/L^. Taking the y-axis to 
lie in the orbital plane, without loss of generality, the Boyer-Lindquist coordinates of the particle at any point of the 
orbit obey the relation 

/(6»,(/>) = cos6' - sine* cose/) tant = (109) 

Symmetry ensures that the self-force has no component perpendicular to the orbital plane and therefore 

pe 

u" ua = F" Ua = =^ = ^ =sin^9Vtan^ L-cot^9 (110) 



28 



In this, Tia = df/dx" is the normal to the orbital plane. The third term in equation (|69|) thus becomes, 

2pgFe^2F^ [tan^ l ~ cot^ 9] pe = 2p^ [tan^ i - cot^ 6] (111) 
But, Lz — Pcj, and = dLz/dr f Appendix [n|l and pU8|) becomes 

^ = 2 (cot^ + tan^ t - cot^ 61) ^ ^2 t&n^ lL^L^ = ^dL^ ,-^^2) 
dr ^ ' dr Lz dr ^ ' 

which brings (|108|) to the form (|102|) . The above equation is valid for both instantaneous and time-averaged fluxes. 
Turning on the spin makes pU8|l depart from its spherical value and consequently l to change. It is easy to deduce 
that some leading order terms (with respect to the spin) will be provided by the third term in H108() . This is because 
the leading order change in Fg is linear in a, as can be deduced by inspecting Ryan's expressions for the self-force 
[ssf. The first term, on the other hand, contributes only at O(a^) and beyond. 

One final approximate Q for generic zoom- whirl orbits has been proposed in Ref. [b^I- As we discussed back in 
Section IlII Dl equatorial zoom- whirl orbits radiate energy and angular momentum at such rates as -E « il^L^- This 
is a manifestation of the fact that a large portion of the orbital period is spent at r « rp, with the small body whirling 
in a quasi-circular path. A similar behaviour is observed for inclined zoom-whirl orbits (in which case ^ Tg,T^), 
hence it is natural to put forward the conjecture that that these orbits also radiate in a 'circular' manner. That 
this could be true is suggested by the Kennefick-Ori formula (|68|l . For a zoom- whirl orbit and for motion near the 
periastron, r « rp, so we should have « 0; consequently, the unknown third term in (|68|l should be negligible. 
Then, the conjecture is that the resulting expression for Q 

Q^n,EE + n,L^L, + 0{u-') (113) 

describes the evolution of the Carter constant for all generic zoom-whirl orbits and with increasing accuracy as the 
orbit approaches the separatrix. We emphasise that this approximation should hold even for orbits deep in the black 
hole's strong-field. This conjecture could become a practical tool once a code that calculates E and Lz for generic 
orbits is developed. Furthermore, a direct comparison with (|102|l should be a useful guide for the accuracy of the 
L = constant rule in strong-field situations. Future computation of the self force will provide the ultimate test for 
both approximations. 



C. Non-equatorial inspirals and waveforms 

Armed with the above approximate expressions for Q it is possible to generate inspirals of generic Kerr orbits. 
These results are particularly interesting, as there are still no analogue Teukolsky-based results. As a warm up we 
first discuss circular-inclined inspirals, for which the circularity theorem fSection IIII E|l ,2^] , ,,32,] , i^Jl guarantees that 
they will evolve as a sequence of circular orbits, and fixes Q. Consequently, for these orbits the rates H87|I .H88 |I become, 

P = H-^A-Lz,cE + E,,Lz) (114) 
^ = H-^c{Lz,pE ~ E^pL,) (115) 
e = e = (116) 

where i/circ = E^^Lz^p — Lz,lE^p. The leading-order expression for p, Z are [4^ . 

64 /X (mV . 244 a ( mV^'^ . ^^^^^ 



Some examples of inclined circular inspirals are shown in Figure lTH and some concrete data are provided in Table HVl 
According to these results, the hybrid approximation outperforms the leading-order calculation. Both approxima- 
tions predict that the inclination angle increases, especially close at the separatrix. However, the increase predicted 
by eqn. H117|) is far too large, particularly for rapidly spinning black holes. The inspiral predicted by the hybrid 
approximation is closer to what is seen in rigorously computed inspirals |l4j . Nonetheless, it too shows an increase in 
i that is probably excessive. After all, holding t constant produces an inspiral sequence that is probably closest of all 
to strong-field calculations and should be acceptably accurate. 

The results presented so far (equatorial, circular-inclined) have established the reliability, as well as limitations of 
the hybrid scheme. Of course, the true purpose of the scheme is the construction of generic inspirals; we discuss them 
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FIG. 14: Circular, inclined inspiral from Ref. |63 |. The black hole spin is a = 0.5M (top) and spin a — 0.9M (bottom). The 
solid lines show inspiral using the hybrid approximation; the dotted lines show the leading order inspiral prediction. The dashed 
curve represents the separatrix. 



next. Generic inspirals are governed by formulae (jHZl , (|HH1) • Before the hybrid approximation, the only available result 
on generic Kerr inspirals was Ryan's weak- field analysis From that work we extract the E,Lz fluxes H9U|) . H91|I . 
Similarly, the natural first choice for Q is eqn. 1)10111 . However, as discussed in Ref. [g^I, the resulting inspirals are 
quite disappointing, particularly in the strong field, unless the orbit is close to the one of the above special cases. The 
root of the problem lies in the use of eqn. H1U1|) . which apparently is not as accurate as we would require it to be. The 
qualitative behaviour of the inspirals is more reasonable when the spherical rule H1U2|I is employed, which also fixes 
Z = 0. The discussion in Appendix 1X1 gives support to the conjecture that this rule should be a good approximation 
for generic inspirals. 

Using the fluxes 19U|I and H91|) together with the Z = rule produces inspirals that agree with the leading-order 
results when p 3> M , smoothly converge to the equatorial results for i ^ and t ^ tt, and that exhibit an e > 
region near the separatrix. All these are features we would expect to be present in the true generic inspirals. Examples 
of generic inspirals are shown in Figs.^l&EI where several interesting features can be seen. The trajectories for 
i < 90° are qualitatively similar to the equatorial, prograde trajectories shown in Fig. 1111 In particular, each such 
trajectory passes through a critical point at which e = after which eccentricity grows. The system typically has 
substantial non-zero eccentricity when it reaches the separatrix. Also, note the excessive growth in eccentricity near 
the separatrix for a = 0.9M and t = 30°. At shallow inclination angle, the separatrix is very deep in the black hole's 
strong field, so the inspiral proceeds to small r before plunging. Just as in the case of equatorial orbits for a = 0.9M, 
the weak-field flux formulae that we use cannot be trusted this deep into the Kerr strong field. 

There are other regions of the {p, e, l} space where the hybrid inspirals exhibit serious pathologies: nearly polar 
orbits {i ~ 90°) and nearly circular orbits (e ~ 0) tend to circularise rapidly, before reaching the separatrix. This is 
not a true physical effect; indeed, it is found that the behavior of inspirals exhibits a rather sharp transition as the 
inclination angle goes from slightly below 90° to slightly above. For nearly circular inspirals we either observe a rapid 
circularisation (for l > 90°) or an unphysically large radius where e changes sign. 

An array of improvements is available for generic inspirals "g?]. First, we can use eqn. H10()|) for I. Second, we 
can build and utilise higher order FN fluxes, instead of l|90|) . H91|l . by combining three different FN results: Ryan's 
fluxes, Tagoshi's fluxes for equatorial orbits with small eccentricity and fluxes for circular-inclined orbits with 
small inclination, derived by Shibata et al |52| . The resulting fluxes can be found in Ref. 

Use of these 'new' 2FN fluxes results in an overall improvement of the inspirals, especially in the strong field regime. 
However, they do not ameliorate the pathological behaviour near t — 90° and e « 0. The root of the problem is found 
when we examine closer the e, p expressions, eqns. (|87|) . For the former we find (a detailed analysis will appear in 
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ajM 


p/M t(deg) 


Calculation 


(M//.)p 


Frac. diff. in p 




Frac. difF. in I 


0.95 


7 62.43 


Numerical 
Hybrid 
Leading order 


-4.657 X 10"^ 
-4.497 X 10"^ 
-2.750 X 10"^ 


0.0344 
0.4095 


1.207 X 10"" 
2.639 X 10"* 
3.080 X 10"'' 


1.1864 
1.5518 


0.5 


10 67.56 


Numerical 
Hybrid 
Leading order 


-1.583 X 10"^ 
-1.645 X 10"^ 
-1.194 X 10"^ 


0.0392 
0.2457 


1.546 X 10"^ 
2.043 X lO"'^ 
2.377 X 10"^ 


0.3215 
0.5375 


0.5 


10 126.76 


Numerical 
Hybrid 
Leading order 


-2.329 X 10"^ 
-2.341 X 10^^ 
-1.414 X 10"^ 


0.0051 
0.3929 


1.892 X lO"'^ 
1.643 X 10"^ 
2.060 X 10"^ 


0.1316 
0.0888 


0.9 


10 74.07 


Numerical 
Hybrid 
Leading order 


-1.544 X 10"^ 
-1.567 X 10"^ 
-1.169 X 10"^ 


0.0149 
0.2429 


2.715 X 10"^ 
3.857 X 10"^ 
4.452 X 10"^ 


0.4206 
0.6398 


0.9 


10 131.57 


Numerical 
Hybrid 
Leading order 


-3.253 X 10~^ 
-3.082 X 10"^ 
-1.545 X 10"^ 


0.0526 
0.5250 


3.887 X 10-5 
2.612 X IQ-'^ 
3.464 X 10-5 


0.3280 
0.1088 


0.5 


6 67.81 


Numerical 
Hybrid 
Leading order 


-2.020 X 10"^ 
-1.779 X 10~^ 
-5.082 X 10"^ 


0.1193 
0.7484 


2.094 X 10-" 
2.992 X lO-"* 
3.954 X 10-"* 


0.4288 
0.8882 


0.9 


6 54.64 


Numerical 
Hybrid 
Leading order 


-7.846 X 10"^ 
-6.950 X 10'^ 
-3.598 X 10^^ 


0.1142 
0.5674 


2.015 X 10-" 
5.486 X 10-" 
6.268 X 10-" 


1.7226 
2.1107 


0.9 


6 99.55 


Numerical 
Hybrid 
Leading order 


-74.32 
-48.02 
-6.593 X 10-2 


0.3539 
0.9991 


6.337 X 10-" 
5.241 X 10-" 
7.580 X 10-" 


0.1729 
0.1961 



TABLE IV: Comparing the rates p, L for Kerr circular-inclined orbits. The accurate numerical data are from p4| . the hybrid 
and leading-order results refer to eqns. (11161 . 1117L respectively. Typically, the hybrid scheme performs much better than the 
leading-order approximation when compared to the rigorous numerical data. The only case which this is not true is for i of 
retrograde orbits. Nevertheless, this small inaccuracy has no impact on the calculation of generic inspirals, as we assume that 
I = constant. Note the enormous difference between the numerical and the leading order results in the Table's final entry. This 
is because that point is fairly close to the separatrix between stable and unstable orbits. Since the leading-order calculation 
has no notion of this separatrix, it is particularly inaccurate here. 
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[7Vi(p, i)E^ 7V2(p, i) + A^ab, Q\ + 0{e) 



The exact fluxes would exactly eliminate the first divergent at e = term, i.e 

7Vi(p, I) E[p, e = 0, i) + N2{p, i) L,{p, e = 0, + mip, i) Q{p, e = 0, = 



(118) 



(119) 



leaving e ~ e which is a well known result for nearly circular orbits [23 , [Ig , Hal ■ It is interesting to note that for 
equatorial orbits, the condition (|119|l reduces to the familiar E = fl^L^ relation. Neither Ryan's original fluxes nor 
the improved 2PN order fluxes satisfy this relation, which explains the bizarre behaviour seen for small eccentricity. 
Looking at H119|) a suitable flx would be. 



^(ncw) ^(old) 



-J-E(p,e^O,i) + -^Q{p,e = 0,L) 



(120) 



while retaining the same E, Q. A similar fix for nearly polar orbits can be devised by studying the limit l 7r/2 of 
e. However, we can arrive at the same result in a more direct way by considering. 
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FIG. 15: Generic inspiral at several inclination angles into a hole with spin a — 0.5M (top pair) and a = 0.9M (bottom pair). 
In all plots, the dashed line represents the separatrix between stable and unstable orbits. As discussed in the text, the inspirals 
are forced to lie in surfaces of constant t. The initial orbital parameters are {pi,ei) — (20M, 0.2), (20M, 0.4), (20M, 0.99), and 
(14M,0.99), and t = 30°, 60°, 80°, 100°, 120°, and 150°. 



Since Q must be well-behaved for t — > 7r/2 it follows that, 



Km 

,^7r/2 



• 2 

sm L 



cos L 



(122) 



This condition is again not exactly satisfied by our Post-Newtonian fluxes, hence they need to be appropriately 
modified. As before, we can just modify as in 1)120(1 . 

Including all available improvements we produce some 'final' generic inspirals, shown in Fie. 1161 These inspirals are 
well behaved throughout the {p, e, t} space, exhibiting a reasonable size for the e > region and for the total change 
in i. They should be fairly good approximation of the true inspirals, but only the comparison with future rigorous 
self force-based inspirals will decide if this statement is true. 

Our discussion on the hybrid approximation concludes with the presentation of a generic kludge waveform (see 
Fig. I17|l . which is calculated according to the prescription given in Section FlV Al The selected orbit is zoom- whirl. 
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FIG. 16: Improved hybrid Kerr generic inspirals, shown on the p-e and p-t planes (left and 
right panels respectively), for a=0.8M spin. Initial orbital parameters are 6^=0.99 and {pi,Li) = 

(20M,60°),(20M,80°),(20M, 100°),(20M, 120°),(14M, 60°),(14M,80°),(14M, 100°),(14M, 120°). Each curve terminates at 
the respective separatrix. These results can be compared with the earlier hybrid inspirals of Fig. 1151 Note the significant 
residual eccentricity and the minor change in inclination. 
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FIG. 17: Kerr generic waveform calculated with the hybrid approximation, using the quadrupole-octupole formula 1)99^ . The 
orbital parameters are p — 3.5M, e = 0.4, t — 40° and the back hole spin a = 0.99M. Note the combined effects of zoom-whirl 
motion and Lense-Thirring precession in the waveform's pattern. 

with p = 3.5M, e = 0.4 and t = 40°, for spin a = 0.99Af. The waveform exhibits a combination of characteristics 
already seen separately in equatorial and circular-incUned waveforms, and is an example of what LISA could actually 
observe. 



V. MAPPING SPACETIMES WITH LISA 

It is universally accepted that the 'dark' objects harboured in galactic nuclei are in fact black holes as predicted 
by General Relativity. But it should also be acknowledged that this belief is partially rooted in our faith to General 
Relativity itself, rather than direct astronomical evidence. So far, it has not been possible to actually probe the 
gravitational field in the vicinity of these dark objects, and determine whether it is described by the Kerr metric or 
something else 's^. Existing astronomical observations can only probe the massive body's asymptotic Keplcrian field 
(there is, however, evidence of existence of event horizons in accordance with the ADAF accretion model 4]). This is 
strong evidence in favor of massive black holes, but is far from being taken as decisive proof of the Kerr metric. 
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LISA will offer a unique opportunity to clarify this matter. Gravitational radiation from EMRI events will serve as 
an optimal probe of black holes spacetimes, verifying their Kerr identity (the most likely outcome of this experiment) 
or identify a more exotic massive object. In order to carry out such a 'spacetime-mapping' program, the LISA data 
analysis community will require the construction of waveform templates for EMRI around a non-Kerr object. Despite 
its potential importance, this program still remains in its infancy. 

The foundations were laid out by Ryan, in a series of papers [iJl > |z3 > Il3 > • The central idea of his analysis is to 
express the general stationary- axisymmetric spacetime in terms of multipole mass and current moments, denoted as 
(£ = 0, 2, 4, ...) and Si {£ = 1, 3, 5, ...), respectively, following the formalism of Geroch, Hansen, Fodor, Hoenselaers 
and Perjes |75| . The metric then can be written in the form, 

ds^ = ~e^+^ dt^ + e^" {dr^ + d9^) + e<-^ sin 6*^(^0 - ujdtf (123) 

The metric functions 7, §j^, ol depend only on r, 9 and each one of them can be written as an expansion in multipole 
moments. For example |74| . 

M 

(5 = V -2^-4^F2«(cos6') + (higher order terms) (124) 

^^^^ fjtA'lX ~\~ X 

n=0 



E 



2 5„_i P2l„_i(cOS0) ^ . . ^ . 

+ (higher order terms) (125) 



2n - 1 r2"+i sin 9 

71—1 

where P2n, ^2^1-1 Legendre and associated Legendre polynomials, respectively. These expansions illustrate the 
lowest order appearance (with respect to 1/r) of each multipole moment: higher moments 'weight' higher powers of 
1/r. The prescription for writing the metric components in terms of multipole moments is given in |7lj |. fj^ and 
makes use of the classic Ernst potential [t^ . 

Kerr spacetime is special in the sense that its higher multipole moments are all 'locked' to the first two moments, 
the mass M = Mq and spin J = Si = aM, 

Ml +iSi ^ M{iaY, £ = 0,1,2,... (126) 

In essence, this statement is the famous 'no-hair' theorem of General Relativity the two-parameter (Af , a) Kerr 
metric is the only possible description for the exterior field of (uncharged) black holes with the corresponding mass 
and spin. 

Given the metric (|123|l . the task is to study geodesic motion of a test-body and the induced gravitational pertur- 
bation. However, neither of these issues is straightforward. The main difficulty that arises is that many of the special 
properties of Kerr spacetime are now lost. Firstly, the fact that the new spacetime (|123|l is not necessarily of algebraic 
type V means that it is not a priori possible to carry out a program that would lead to a Teukolsky-like decoupled 
perturbation equation. Recall that decoupling of the various perturbation equations for the Weyl scalars (in the 
context of the Newman-Penrose formalism) is a property crucially related to the type V character of a spacetime 15^ . 
Secondly, even at the simple level of geodesic motion, the existence of a 'third' constant, like Q, is not guaranteed 
anymore. This is a familiar situation in the study of orbital motion in galactic potentials |77j and it could well lead 
to undesired effects like chaotic behaviour. 

These issues of principle forced Ryan to consider only special cases: (i) the gravitational wave-phase for a circular 
equatorial orbit, and (ii) the orbital frequencies for circular orbits with small eccentricity and small inclination. For 
the former calculation, Ryan resorted to Thome's multipolar-expansion for the flux '78l and obtained. 
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(127) 



where v = ^ M /p. Then the wave-phase can be written as. 



$(t) = 2tt j dt f{t) ^2tt j df f with / = 2 n^/2Ti = n^/n (128) 

For studying geodesic motion one can again write a Hamilton- Jacoby equation which, however, is only separable with 
respect to t and coordinates (a consequence of the spacetime's stationarity and axisymmetry) ; separability with 
respect to r and 9 is no longer possible (apart from the trivial case 9 ~ tt/2). Then, the small body's momentum is 
given by: 

p' = -g't E + g"^ L,, p* = g^^L,-Eg"^ (129) 
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and 



(130) 



Following Ryan, assuming a small eccentricity - small inclination orbit allows eqn. (|130|l to decouple (with respect to 
r and 9) at leading order in e, i. In this way, it is possible to work out the fundamental orbital frequencies jTlj . 



and 



^0 — 



VLm dE 



1 - 3u^ 



2 Si 



45i 
M2 



3M; 



M2 



2M3 



2 4 



7^2 



3 Ma 



1 4 20^1 

2 " 9 M2 



M4 M3 



27 Af2 



+ O(w^) 



(131) 
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If either the phase H128() , H133() or the frequencies (|131|) , H132|l are measured, then it will be possible for LISA to extract 
the first few moments with sufficient accuracy, as suggested by Ryan's analysis jOljIZ^- However, much more detailed 
work is required in order to construct actual templates for waveforms generated by bodies in generic orbits of the 
metric (|123|l . 

The multipole moment formalism, despite its general and elegant character, has several serious disadvantages. 
Expanding the metric as in eqn. H125(l is in practise an expansion in 1/r. For any practical implementation this 
expansion has to be truncated. The induced error is not significant at large values of r but in the strong field regime 
several moments need to be included. 

Serious difficulties emerge from the non-separability of the Hamilton-Jacobi equation and the lack of a 'Teukolsky- 
like' perturbation formalism. The former somewhat complicates the description of geodesic motion, while the latter 
seriously inhibits any attempt to calculate gravitational waveforms, fluxes and back reaction to the orbit. PN results, 
as the ones given by Ryan, are unfortunately not accurate enough for constructing waveform templates for LISA. 

These problems provide sufficient motivation for approaching the problem in a different manner. There is some very 
recent work along this line [t^jUSI which is based on the idea of a 'quasi-Kerr' metric: the underlying assumption is 
that the spacetime is not exactly described by the Kerr metric, but rather by a metric which deviates only slightly 
from Kerr. 

The actual prescriptions for building a quasi-black hole metric is markedly different in the existing studies. In 
Ref. [t^ , Collins & Hughes construct (as a first stab to the problem) a 'bumpy' Schwarzschild black hole by artificially 
adding a certain distribution of mass fl <^ M around a Schwarzschild black hole, in the form of (i) point masses at 
the hole's poles or (ii) an equatorial ring encircling the black hole. In this way, the perturbed metric is axisymmetric 
and stationary and has a non-zero quadrupole moment (as well as higher mass moments) , whereas the Schwarzschild 
spacetime has only a mass moment. Collins & Hughes write the resulting metric in the following Weyl form (using 
cylindrical coordinates {pjZ, (/)}), 



(134) 



with '0 = ■00 + "01 7 7 ■ 
induced perturbation. 



= 7o + 7i- The values ■0Oi7o correspond to the exact Schwarzschild metric while 0'i,7i is the 
For example, for ipi we have, 
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The point masses are located at z ± 6 and the ring at p = b. Knowledge of ip immediately provides 7, since they are 
related by the vacuum Einstein equations (see [t^ for relevant formulae). Both 0i, 71 scale as ^ i^/M, and are to be 
kept only to leading order in fi. This operation is equivalent to truncation of the multipolar structure of the metric 
(|134|l to the level where only the mass and quadrupole moment are taken into account. 

Collins & Hughes also worked out the equatorial geodesic equations of motion for the metric 113411 . They find. 
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FIG. 18: Periastron precession angle /S.<j) due to a combined quadrupolar distortion, eqns. I|185|l . I|136|l of Schwarzschild space- 
time. The top curve refer to strong-field results while the bottom curve is the weak-field Keplerian estimate 11381 . 



An orbital parameter which is clearly imprinted on the waveform is the periastron advance A0. Due to the non-zero 
quadrupole moment, geodesic motion in the metric H134|) will exhibit an additional precessional motion in addition 
to the familiar Schwarzschild precession. It is easily calculated from eqns. (|137|l and is plotted in Figure E| (for 
/2 = 0.01Af andfe = 2.5M). 

In the weak-field limit (p ^ M) the periastron precession takes the simple form. 
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According to H138|) , a prolate (oblate) quadrupole deformation will induce a positive (negative) precession in agreement 
with Newtonian expectations. As Fig. ^| suggests the weak-field result seriously underestimates At/) for p ^ 20M. 
Either way, the perturbation ?Ai,7i can impart a significant correction to the small-body's orbit. 

In their study 80], Glampedakis & Babak build a quasi-Kerr spacetime in a somewhat more natural way, by 
combining the exact Kerr metric with the celebrated Hartle-Thorne (HT) exterior metric [sj. This latter metric 
describes the spacetime outside any arbitrary, but slowly rotating, stationary and axisymmetric body. It is accurate 
up to 0{J^) where J is the body's spin. The HT metric depends on three parameters: the mass M, spin J and 
quadrupole moment Q. For the special case, 
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it reduces to the ©(a^) Kerr metric. We can then assume an approximation where the quadrupole moment has a 
slightly 'wrong' value. 
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kcrr 



with e ^ 1. Then the Hartlc-Thornc metric can be written. 
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The quadrupole term e hab represents the leading order deviation from the Kerr metric. The ansatz adopted in |80i | is 

to replace the 0{a^) Kerr metric with the full Kerr metric, i.e g^^ sib^'^"''- However, prior to this operation, 

one more manipulation is required in (|141|l . The Kerr limit (e 0) of the original HT metric (H^) does not yield 
the Kerr metric in its Boyer-Lindquist coordinate form ^^l- This not surprising, as the original set of coordinates 
used by Hartle & Thorne are not Boyer-Lindquist. On the other hand, there is a strong reason for which we would 
like this limit to be satisfied, and is related to the separability of the Hamilton-Jacobi and wave equation. The Kerr 
limit of the HT metric leads to non-separable dynamical equations (such as the Hamilton-Jacobi and the scalar wave 
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equation). These equations are separable provided that we are deahng with a type V spacetime (which is true for 
the HT metric only when e ^ ) and for a certain class of coordinate frames |23|,|2J|. As the Boyer-Lindquist 
frame is included in this privileged class, it makes sense to express the HT metric in these coordinates. In fact, the 
relevant transformation from Boyer-Lindquist to Hartle-Thorne coordinates was given in the original paper |8lj| and 
it is trivially inverted, at 0{J^) accuracy 80]. 

Putting all the pieces together, the ansatz for building a quasi-Kerr metric in Boyer-Lindquist coordinates is: 

ffa6=5if"^+e/ia6 (142) 

where, 

/i" = {\-2Mlr)-\{l-?,cof?e)Ti{r)l h" = {I - 2M lr)[{l ~ icof? 9) Ti{r)] (143) 

h'' = -\[{l-icos'e)T2{r)l h-^-P = —-}—[{l-icos'e)T2{r)] (144) 

r'^ sm 9 

h"^ = (145) 

The radial functions J^i,2{i") are written explicitly in Appendix IfI 

Given the metric H142|) . it is straightforward to obtain equations of motion for equatorial orbits, 

P' = pW'-„)+'TT7W R2 = -r(.--2M)|MV^, + i?(^,-^,)] 



and Ro ^r^E^ ^ r{r ~ 2M) [fi'^r^ + hi] (147) 



where P^kcit) pure Kerr four- momentum (eqns. Q). We have to point out that if we were merely interested in 

equatorial motion only, then shifting from Hartle-Thorne to Boyer-Lindquist coordinates is an unnecessary step, as 
the Hamilton Jacobi equation is trivially separable when 9 = n 12. However, the coordinate transformation is required 
in the general case in order to ensure that the Hamilton-Jacobi equation is separable in the e ^ limit. Then, one 
can still obtain approximate equations of motion for generic orbits by means of canonical perturbation techniques 
(see ^3 for treatment of the perturbed Kepler problem); these equations will appear in Ref. (80| . 

The main goal of the quasi-Kerr framework is to facilitate the computation of strong-field fluxes and waveforms. 
So far, the only relevant results have been obtained by applying the hybrid approximation for the metric H142I) . in 
combination with H148|l |80l | (indeed, it would be straightforward to repeat the same calculation for the metrics (|123|l . 
(|134|l . always assuming equatorial motion). Some preliminary results are shown in Figs. 1191 fc 1201 where we compare 
quasi-Kerr versus exact Kerr waveforms and quantify the difference of the periastron shift. These results suggest that 
(say) a 10% deviation in the quadrupole moment would lead to significant changes in the orbit and waveform. 
Possible non-Kerr character of black hole spacetimes is a factor that should be taken into account in building data 
analysis tools for LISA'^. 

The computation of rigorous fluxes, waveforms and backreaction to the orbit is the major challenge that lies ahead 
for the quasi-Kerr program. The metric H142() has the clear advantage of being arbitrarily close to type T). The hope 
is that this property may eventually allow the formulation of a 'quasi- Teukolsky' equation for the perturbed Weyl 
scalars, in which case, one will be able to obtain results in a fashion similar to the one described in Section IIIII In the 
unfortunate case where this scenario will not flourish, we may be forced to deal with direct metric perturbations. These 
are described by ten coupled equations (plus four gauge constraints) which will require solution in the time-domain 
(as a '2-1-1' problem). Clearly, this is a much more complicated task compared to the solution (either in the time or 
frequency domain) of the single Teukolsky equation. It is interesting to note that recent approaches to the calculation 
of the gravitational self-force, is via time-evolution of black hole metric perturbations (see other contributions in the 



In this context, an important issue that needs to be addressed is whether a Kerr waveform with a different set of parameters could 
match the quasi-Kerr waveform. 
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FIG. 19: Quasi-Kerr (solid curves) versus Kerr waveforms for the orbit p = lOM, e = 0.5, t = 0°. The black hole spin is 
a = 0.5M and the Kerr deviation parameter e = 0.1. In calculating the Kerr waveforms we have either ignored or included 
the radiative orbital evolution (setting ^/M — 10~^ ), dashed and dashed-dotted curves respectively. To demonstrate how 
promptly the quasi-Kerr waveform dephases we show two different time windows: an early adiabatic stage and a much later 
{t ~ Trr) stage where backreaction effects have become noticeable. Note that the quasi-Kerr waveform already deviates from 
its Kerr counterpart after only few orbits, where the two Kerr waveforms are essentially indistinguishable. 
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FIG. 20: Number of orbits required to accumulate Ai^ = 7r/2 periastron advance for equatorial orbits in the quasi-Kerr metric 
11421 1 as a function of e. The spin is a = 0.5M and the orbital paramaters are: p = lOM, e = 0.5 (solid curve), p = 15M, e = 0.5 
(dashed curve). 



present Volume). Existing calculations have been focused on the Schwarzschild case; generalisation to Kerr is at 
early stages. It is likely that the experience gained by these calculations could be extremely useful for the quasi-Kerr 
problem. 



VI. RECOILING BLACK HOLES 



Gravitational radiation generated by an inspiralling and subsequently coalescing binary system is not only a carrier 
of energy and angular momentum but also of linear momentum. As a result, the binary's center of mass experiences a 
'recoil' kick [69^ . This recoil effect is expected to be of potential astrophysical importance in situations where massive 
black holes are involved. At the top of the scale we have the inspiral and final merger of supermassive black holes, 
as the aftermath of galactic mergers. At a smaller scale, there may be intermediate mass black holes merging at the 
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centers of globular clusters. 

In particular, it is important to obtain an estimate for the final net recoil velocity of a merged black hole binary, and 
compare it to the escape velocity of the host stellar system. In this problem one has to consider a binary composed of 
comparable mass black holes. The general formula for the averaged rate at which linear momentum is removed from 
the system (which is the opposite of the amount carried by gravitational waves) is [2Q| |. 



pk ^ _ f 
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16^ 



(149) 



where n*^ is the radial unit vector. Unlike the energy and angular momentum fluxes which first appear at the level 
of the source's mass quadrupole moment, linear momentum is a higher order effect and requires combination of 
quadrupole mass and current moments together with mass octupole moment. This leading order formula was first 
derived by Bekenstein ,63|| , 
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where I*^ , Al*^*^, 5*^ the symmetric and trace-free mass parts of 7*^ , S'*^'', M*^*^, see Ref. ^3 fo'" details. 

A first recoil result was obtained by Fitchett [s^, using H15U|) under the assumption of Kcplcrian motion for the 
binary. For circular orbits Fitchett 's result is, 
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where q = M1/M2 < 1 is the binary's mass ratio, rtcim is a 'cut-off' radius where gravitational wave emission 
terminates, and f{q) = q'^{l ~ q^)/{l +q)^- According to H151|) . for the expected rterm ~ few x M, the recoil velocity 
can be vp ^ 10—100 km/s. Typical escape velocities span a range ~ 10 — 10"^ km/s, from globular clusters to galaxies. 
In view of these numbers, Fitchett's result clearly establishes the astrophysical importance of the recoil effect and 
provides motivation for more detailed work. 

More detailed work means taking into account strong-field effects in the motion of the coalescing black holes since, 
after all, the momentum flux grows with decreasing separation for the two black holes. However, this is a venture into 
the thorny area of fully nonlinear General Relativity hence more suitable for hard-core numerical relativity (see [s^l for 
estimates on binary black hole recoil) . Fortunately, some valuable input can be provided by much simpler perturbative 
methods. Although black hole perturbation theory (such as the Teukolsky formalism) assumes an extreme mass ratio 
g <C 1, results can be extrapolated to (say) q ^ 0.5 with a modest error (with the same reasoning we can also neglect 
the spin of the lighter black hole). In fact, the recent study by Volonteri et al [s^l suggests that low-redshift ( z < 1) 
black hole mergers predominantly occur with g ^ 0.1. 

Within the Teukolsky framework the first relevant calculations were performed by Fitchett & Detweiler [s^ for 
circular Schwarzschild orbits. Linear momentum fluxes were also computed for bound and plunging parabolic orbits 
in Kerr and Schwarzschild spacetimes |40l |. 

In terms of the amplitudes appearing in the expansion for ?/'4, 
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the linear momentum flux is 
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Note that these expressions are valid only for circular equatorial orbits. The total recoil velocity simply follows by 
integrating (|155|l from some initial time to (when the binary separation is large and ^^.^(to) = 0); 
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FIG. 21: Recoil velocity as a function of black hole spin a/M for q ~ 0.1. (a) Recoil velocity up to the separatrix. The 
sohd curve is the Teukolsky-based result while the dashed curve shows the Newtonian recoil prediction eqn. 11511 . which is 
substantially higher for a > 0.5M prograde inspirals. (b) Upper and lower limits for the total recoil. The shaded region 
represents the uncertainty in the final kick velocity. The upper limit curve is generated by assuming a continuation of the 
power-law behaviour Ptot ~ p~" that fits the inspiral down to the separatrix. On the other hand, the lower-limit curve is 
generated with the help of the hybrid scheme. 



with Pa; J, calculated for a sequence of quasi-circular orbits. The most recent study on recoilling black holes is the one 
by Favata et al |86j who focused on Kerr circular equatorial orbits. In a realistic scenario there is no a priori reason 
to expect the orbit to be 'equatorial' (for a moderate mass ratio system, the term equatorial should be interpreted as 
alignment of the total orbital angular momentum and the spin of the heavier black hole) . Still, one can gain useful 
information by examining the two extreme opposite cases of equatorial prograde and retrograde motion. 

Estimating the recoil from eqns. (|155|l . (115611 is an accurate procedure (apart from the error introduced by the finite 
mass ratio) as long as the orbit is stable, p > ps- Results on the recoil velocity, as computed by Favata et al, are 
shown in Figure f?n According to this data, Fitchett's leading order formula is in fairly good agreement with the 
Teukolsky-based results. Some noticable deviation appears for prograde motion around a > 0.5 A/ Kerr holes. The 
orbit ventures deep inside the hole's strong field rendering the leading-order result (|151|l unreliable. 

The major uncertainty in the recoil calculation comes from the plunging part of the inspiral. Indeed, this part 
could well provide the dominant contribution to the final recoil velocity, especially when the motion is retrograde (in 
which case the separatrix is beyond 6M). In contrast, prograde orbits with a ~ M receive little contribution from 
the plunge since ps ~ r^. 

Somewhat surprisingly, there is no available Teukolsky-based analysis for this kind of orbits. In Ref. ^8(^ two 
different approximations are employed in order to obtain some rough estimate for the contribution of the plunge 
portion of the inspiral. Firstly, a power-law behaviour P ^ is assumed, as a continuation of the similar behaviour 
during the adiabatic quasi-circular inspiral (see for example Fig. taken from Ref. 'sS] ) . The resulting recoil velocity 
is designated as 'upper limit' in Fig. 1211 Secondly, the hybrid approximation fSection llVp can be employed by using 
the exact plunging geodesic trajectory in the leading-order formula H15U|) . In this way, the 'lower-limit' recoil curve 
of Fig. 1211 is generated. The hybrid result should not be such a bad approximation given the overall accuracy of 
this calculation. This is what Fig. [53 suggests, where we compare Phybrid against Ppeuk for p > pg. In any case, 
what we learn with confidence from Fig. 1211 is that for a merging binary with the heavier member rapidly spinning 
and the lighter one in prograde motion, the resulting recoil velocity is ~ 50 — 100 km/s smaller than the Newtonian 
expectation. 
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FIG. 22: Linear momentum flux at infinity for Kerr circular-equatorial orbits. We compare results obtained by the hybrid 
approximation (solid and dashed curves for prograde and retrograde motion respectively) with accurate Teukolsky-based results 
(box and circle symbols) and the leading order prediction (dotted curve) as given by eqn. 115111 . The black hole spin is set to 
a = 0.9M. 

VII. CONCLUDING REMARKS 

We have presented a variety of material on EMRIs, from rigorous Teukolsky-based results (both in the frequency 
and time domains) for special orbits (equatorial-eccentric and circular-inclined) to approximate results for generic 
orbits (hybrid approximation). Guided by these studies we have a fairly detailed knowledge of the gravitational 
dynamics of the system and of the associated waveform. However, there is still plenty of distance until the final goal 
of developing realistic waveforms for LISA is reached. Presumably, this program will be fulfilled in the near future 
by the self- force framework. But even within the limits of the standard flux-balance/Teukolsky formalism there is 
much work to do. A first priority is the computation of rigorous waveforms and fluxes from a test-body in a generic 
orbit, without worrying about the orbital evolution and Q. Fortunately, this analysis is expected to appear soon 
(3 lj| . Another important task is the construction of Teukolsky-based inspirals, as the ones discussed in Section IlII Fl 
We already have the tools for producing equatorial-eccentric inspirals, and it will even be possible to 'fake' generic 
inspirals by using rigorous fluxes and some approximation for Q (as the ones described in Section IIVB|I . 

Another layer of complexity can be added regarding the modelling of the small body itself. The assumed 'point- 
particle' is in reality another Kerr black hole (or some other compact body) with intrinsic properties beyond the mass 
/i. Modelling the small body, one could take one more step and include the spin. Formal equations of motion for 
a spinning test-body have been known for a long time "s^, but only recently there were actual calculations for the 
particular case of Kerr spacetime |9nj |. The coupling of the body's spin with the background curvature results in 
a non-geodesic motion with extra spin degrees of freedom and loss of one integral of motion (the Carter constant). 
Additional work is required in order to understand the motion of a spinning body under the influence of radiation 
reaction (plus the backreaction to its spin), within both the adiabatic flux-balance and the self- force framework. 

If we are prepared to accommodate the idea that LISA might actually discover non-Kcrr massive objects, then we 
must be prepared to study EMRIs around non-Kerr metrics. The material of Section is just a first modest step 
towards that direction. The actual computation of 'realistic' waveforms/orbital evolution could develop to a daunting 
task, likely to involve direct time-evolution of metric perturbations. Also it is not clear yet if by including only the 
leading-order quadrupole moment deviation is sufficient for actually telling the difference between a true quasi-Kerr 
object and a Kerr black hole with slightly different set of M, a. If this turns out to be the case, and additional 
multipole moments are required, we may then have to revitalise (and perhaps modify) Ryan's multipolar formalism. 
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APPENDIX A 



The occurrence of a 'third' orbital constant Q in axisymmetric gravitational fields is not an exclusive feature of 
general relativity. For example, it is familiar from Newtonian celestial mechanics applied to orbital motion in galactic 
gravitational potentials (see, for example (77], where the third constant is denoted /). The departure of Q from 
+ can then be attributed to the 'asphericity' of the potential. If such a potential does not deviate very much 
from sphericity, (the square of the total angular momentum) turns out to be almost constant, so that Q should 
be, after all, nearly I? — L\. 

It is straightforward to check whether this behaviour of occurs in Kerr spacetime. The definition we use for 1? 
is identical to that used in Schwarzschild spacetime, 

i2=p2 + (sin0)-2p2 ^ (Al) 

For the Carter constant Q we have. 
Combining these two expressions gives 

Q = - LI + cos^ 9{f? - E^) . (A3) 

In other words, Q can be interpreted as the projection of the total angular momentum on the equatorial plane, apart 
from the 'aspherical' term cos^ 0(/x^ — E"^). This interpretation makes sense when the aspherical term is small, 
that is, when a ^ M (slow rotation) and/or E ^ 1 (weak-field orbits). In practice, we find that this term is often 
significantly smaller than the preceding terms even for motion in strong-field regions of rapidly rotating holes. We 
illustrate this in Fig. |2S1 showing how the quantity SL^ = L"^ / [Q + i^) — 1 varies with time for a variety of generic 
orbits around a rapidly spinning hole. 

Examining Fig. [531 we see that deviates very little from Q + even when the small body is deep in the 
black hole's strong field; in this sample, the difference is no more than about 1%. This shows that interpreting Q 
as a squared projection of angular momentum into the equatorial plane is sensible. Because Q + is a constant 
quantity, this figure also demonstrates that is nearly constant. This is exactly what we expect for motion in an 
axisymmetric potential that is almost spherical. These pieces of evidence suggest that the Kerr spacetime is not as 
aspherical as we might have expected, at least for the purposes of this argument, lending credence to our suggestion 
that the 't = constant' assumption should be reliable. 



APPENDIX B 

Here we list the integrals appearing in the expressions for the Kerr fundamental frequencies, eqns. (|18(I . H19(I . We 
have 

A = K{k) y(rp, ra) + a" z\ \K(k) - E{k)] X{rp, r^) (Bl) 



with 



and 



Kik) = / , . ■> ^ E{k)^ d^A/l-fcsin^V (B2) 

Jo \/l — k sin ip Jq 

U{zl,k) ^ (B3) 

Jo (1 - sin^ V) V 1 - sin ?A 
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FIG. 23: Examining our notion of "total angular momentum" for strong- field Kerr black hole orbits. Each panel compares the 
angular momentum squared = pg + (sin&)^^p| to Q + L^: the quantity plotted is 5L^ = /{Q + Lj.) — 1. The top panel 
shows these quantities over an orbit with p = lOM, the center panel the quantities over an orbit with p — 5M, and the bottom 
over an orbit p — AM. In all cases, the orbits have eccentricity e = 0.7, inclination l = 60°, and are about a hole with spin 
a = 0.999M. Even deep in the strong field, differs very little from Q + — the greatest deviation in this sample is about 
1%. Since Q + is a constant by definition, this also shows that is approximately conserved over the orbit. 



APPENDIX C 



In this Appendix we provide a brief discussion on the Kennefick-Ori formula for the Carter constant flux Q, see 
Ref. The starting point is to express all three constants C = {E,Lz,Q} in a functional form C — C{x'',pa)- In 
particular, for E, Lz we simply have. 



E ^ -pt, 

while for Q we have a choice between the two (equivalent) forms. 



Q 



cos^ 



sm 



giO,E,Lz) 



Q = [ptir^ + a^) + ap^] - (p^ + aptf - - Apj = W(r, E, L,) - Apj 

Defining the self-force as Fa = Dpa/Dr, Kennefick & Ori derive, 

dC _ dC ^ 
dr dpa 



(CI) 

(C2) 
(C3) 

(C4) 



for the instantaneous rates. Clearly, Ft — ~dE / dr and F^ = dLz/dr. For Q we obtain the twin expressions, 

Q ^ nE^+ n.L^ ^ - 2 A PrFr (C5) 

Q ^ g.E^ + g,L.^+2pgFe (C6) 
dr dr 
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APPENDIX D 



This appendix contains explicit expressions (in terms of E, Lz, Q and their derivatives) for the various functions 
appearing in the formulae H87|) for the rates p, e, I. First, 

H = Q^p £',e -^2,1 ~ Q,p £',t ^2,e ^ Q,e E^p L^^l + Q ,e -£',t -^2,p + Q ,l E^p L^^e ^ Q ,l -£',e ^z,p (Dl) 

and 

bp — Q,LLz,e ~ Q,eLz,n — Lz^iQ^p — Q.^Lz^p (D2) 



— E ^Q fy E Cq — QjiE p E 



dp — E^f.Lz^u ~ E^^Lz^e, — E^^Lz^p ~ E^pLz^u 



(D3) 
(D4) 



APPENDIX E 



The functions /i(e) — /8(e) and .gi(e) — 59(e) appearing in the flux- formulae (|90() . 191|) . H95|I . (|96|I are: 
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APPENDIX F 



The radial functions J-'i.2ir) appearing in the quasi-Kerr metric (|142(l are: 

5(r--M) . , ^^^^ 15r(r-2Af) 



Ti(r) = ^ ^ (2^^^ + 6A/r - 3r^) i ^ — ^ In (Fl) 

8Afr(r-2M)^ ^ 16M2 \r - 2M J ^ ' 

T^(r) = (2Af2 - 3Mr - ir^) + — ?^(r2 - 2A'/2) I ^ ) (F2) 



[1] Danzman K et al, LISA - Laser Interferometer Space Antenna, Pre-Phase A Report, Max-Planck-Institute fur Quantenop- 
tic, Report MPQ 233 (1998) 

[2] Material on a variety of LISA and EMRI topics can be found at the following WWW address: 

http:/ / www. tapir, caltech. edu/listwgl / 
[3] Rees M J in Black Holes and Relativistic Stars, (ed.) Wald R M (Chicago: Univ. Chicago Press) (1998) 



44 



[4] Ichimam S Astroph. J. 214, 840 (1977); 

Narayan R and Yi I Astroph. J. 428, L13 (1994); 

Abramowicz M A, Chen X, Kato S, Lasota J-P, Regev O Astroph. J. 438, L37 (1995); 
Narayan R and Yi I Astroph. J. 452, 710 (1995); 

Narayan R, Garcia M R and McClintock J E Astroph. J 478, L79 (1997) 
[5] Rees M J Class. Quantum Grav. 14, 1411 (1997); 

Sigurdsson S and Rees M J Mon. Not. R. Astron. Soc. 284, 318 (1997); 

Sigurdsson S Class. Quantum Grav. 14, 1425 (1997) 
[6] Freitag M and Benz W Astron. & Astroph. 375, 711 (2001); 

Freitag M Class. Quantum Grav. 18, 4033 (2001); 

Freitag M Class. Quantum Grav. 583, L21(2003); 

Gair J R Barack L Creighton T Cutler C Larson S L Phinney E S and Vallisneri M Class. Quantum Grav. 21, S1595 

(2004); 

Hopman C and Alexander T, preprint, astro-ph/0503672] 
[7] Rees M J Rev. Mod. Astron. 2, 1 (1989) 
[8] Chakrabarti S K Phys. Rev. D 53 2901, (1996); 

Chakrabarti S K Astroph. J. 411, 610 (1993) 
[9] Narayan R Astrophys. J. 536, 663 (2000) 
[10] Chandrasekhar S, The Mathematical Theory of Black Holes, Oxford Univ. Press (1983) 
[11] Cutler C, Kennefick D and Poisson E Phys. Rev. D 50, 3816 (1994) 
[12] Bardeen J M, Press W H and Teukolsky S A Astroph. J. 178, 347 (1972) 
[13] Glampedakis K and Kennefick D, Phys. Rev. D 66, 044002 (2002) 
[14] Hughes S A Phys. Rev. D 61, 084004 (2000); 

Hughes S A Phys. Rev. D 63, 064016 (2001) 
[15] Schmidt W Class. Quantum Grav. 19, 2743 (2002) 
[16] BombeUi L and Calzetta E Class. Quantum Grav. 9, 2573 (1992) 
[17] Ori A and Thorne K S Phys. Rev. D 62, 124022 (2000) 
[18] O'Shaughnessy R Phys. Rev. D 67, 044004 (2003) 
[19] Goldstein H Classical Mechanics (Addison- Wesley, 1950) 

[20] Misner C W Thorne K S and Wheeler J A, Gravitation (Freeman, San Francisco, 1973) 
[21] Carter B Commun. math. Phys. 10, 280 (1968) 

[22] Wald R M General Relativity (The University of Chicago Press, 1984) 
[23] Teukolsky S A Astroph. J. 185, 635 (1973) 

[24] Stewart J M and Walker M Proc R. Soc. Lond. A 341, 49 (1974) 

[25] Mino Y, Sasaki M, Shibata M, Tagoshi H and Tanaka T Prog. Theor. Phys. Suppl. 128, chapter 1 (1998) 

[26] Landau L D and Lifshitz E M The classical theory of fields (Addison- Wesley and Pergamon, 1971) 

[27] Teukolsky S A and Press W H Astroph. J. 193, 443 (1974) 

[28] Hawking S W and Hartle J B Gommun. Math. Phys. 27, 293 (1972) 

[29] Drasco S and Hughes S A Phys. Rev D 69, 044015 (2004) 

[30] Mino Y Phys Rev D 67, 084027 (2003) 

[31] Hughes S A, Drasco S, Flanagan E E and Franklin J, preprint gr-qc/050401^ 

Drasco S and Hughes S A, in preparation (2005) 
[32] Detweiler S L Astroph. J. 225, 687 (1978) 

[33] Sasaki M and Nakamura T Prog. Theor. Phys. 67, 1788 (1982); 

Sasaki M and Nakamura T Phys. Lett. 89A, 68 (1982) 
[34] Breuer R A, Ruffini R, Tiomno J and Vishveshwara C V Phys. Rev. D 7, 1002 (1973); 

Breuer R A and Vishveshwara C V Phys. Rev. D 7, 1008 (1973); 

Chrzanowski P L and Misner C W Phys. Rev. D 10, 1701 (1974) 
[35] Apostolatos T, Kennefick D, Ori A and Poisson E Phys. Rev. D 47, 5376 (1993) 
[36] Kennefick D Phys. Rev. D 58, 064012 (1998) 
[37] Peters P C and Mathews J, Phys. Rev. 131, 435 (1963) 
[38] Ryan F D Phys. Rev D 52, R3159 (1995) 
[39] Detweiler S L and Szedenits E Astroph. J. 231, 211 (1979) 
[40] Oohara K and Nakamura T Phys. Letters 94 A, 349 (1983); 

Oohara K and Nakamura T Prog. Theor. Phys. 71, 91 (1984); 

Kojima Y and Nakamura T Prog. Theor. Phys. 71, 79 (1984); 

Kojima Y and Nakamura T Prog. Theor. Phys. 72, 494 (1984) 
[41] Tanaka T, Shibata M, Sasaki M, Tagoshi H and Nakamura T Prog. Theor. Phys. 90, 65 (1993) 
[42] Shibata M Phys. Rev. D 50, 6297 (1994) 
[43] Finn L S and Thorne K S Phys. Rev. D 62, 124021 (2000) 
[44] Shibata M Prog. Theor. Phys. 90, 595 (1993) 
[45] Kennefick D and Ori A Phys. Rev. D 53, 4319 (1996) 
[46] Ryan F D Phys. Rev D 53, 3064 (1995) 
[47] Hughes S A Phys. Rev. D 64, 064004 (2001) 



45 



[48] Wilkins D C Phys. Rev. D 5, 814 (1972) 

[49] Bckcnstcin J D and Schiffer M Phys Rev D 58, 064014 (1998) 
[50] Press W H and Tcukolsky S A Nature 238, 211 (1972) 
[51] Hartle J B Phys Rev D 8, 1010 (1973); 

Hartle J B Phys. Rev D 9, 2749 (1974) 
[52] Shibata M, Sasaki M, Tagoshi H and Tanaka T Phys. Rev. D 51 , 1646 (1995) 
[53] Poisson E Phys. Rev. D 55 , 639 (1997) 

[54] Davis M, Ruffini R, Press W H and Price R H Phys. Rev. Lett. 27, 1466 (1971) 
[55] Lopez-Aleman R PhD Thesis, Perm State University (2001) 

[56] Lopez-Aleman R, Khanna G and PuUin J Class. Quantum Grav. 20, 3259 (2003) 
[57] Khanna G Phys. Rev D 69, 024016 (2004) 
[58] Martel K Phys. Rev. D 69 044025 (2004) 

[59] Krivan W, Laguna P, Papadopoulos P and Andersson N, Phys. Rev. D 56 3395 (1997) 

[60] Martel K and Poisson E Phys. Rev. D 66, 084001 (2002) 

[61] Campanelli M and Lousto C O Phys. Rev. D 59, 124022 (1999) 

[62] Lousto C O and Price R H 1997 Phys. Rev. D 56 6439 

[63] Laguna P, talk given in the 7th Capra Radiation Reaction Meeting (2004) 

[64] Glampedakis K, Hughes S A and Kennefick D Phys. Rev. D 66, 064005 (2002) 

[65] Sasaki M Prog. Theor. Phys. 69, 815 (1983); 

Ruffini R and Sasaki M Prog. Theor. Phys. 66, 1627 (1981) 
[66] Tagoshi H Prog. Theor. Phys. 93, 307 (1995) 
[67] Gair J R and Glampedakis K, in preparation (2005) 

[68] Babak S V, Fang H, Gair J R, Glampedakis K and Hughes S A, in preparation (2005) 
[69] Bekenstein J D Astroph. J. 183, 657 (1973) 

[70] Press W H, Phys. Rev. D 15, 965 (1977) 
[71] Ryan F D Phys. Rev D 52, 5707 (1995) 
[72] Ryan F D Phys. Rev. D 56, 1845 (1997) 
[73] Ryan F D Phys. Rev. D 56, 7732 (1997) 
[74] Ryan F D Phys. Rev D, 55, 6081 (1997) 
[75] Geroch R J. Math. Phys. 11, 2580 (1970); 
Hansen R O J. Math. Phys. 15, 46 (1974); 

Fodor G, Hocnsolaers C and Porjes Z J. Math. Phys. 30, 2252 (1989) 

[76] Ernst F J Phys. Rev, D 167, 1175 (1968) 

[77] Binney J and Tremaine S Galactic dynamics (Princeton University Press, Princeton, New Jersey, 1987), chapter 3. 

[78] Thorne K S Rev. Mod. Phys. 52, 299 (1980) 

[79] Collins N A and Hughes S A Phys. Rev. D 69, 124022 (2004) 

[80] Glampedakis K and Babak S V, in preparation (2005) 

[81] Hartle J B and Thorne K S Astroph. J. 153, 807 (1968) 

[82] Fitchctt M J Mon. Not. R. Astron. 203, 1049 (1983) 

[83] Andradc Z and Price R H Phys. Rev D 56, 6336 (1997); 

Anninos P and Brandt S Phys. Rev. Lett. 81, 508 (1998); 

Brandt S and Anninos P Phys. Rev. D 60, 084005 (1999); 

Lousto C O and Price R H Phys. Rev. D 69, 087503 (2004) 
[84] Volontcri M, Madau P, Quatacrt E and Recs M J Astroph. J. 620, 69 (2005) 
[85] Fitciict M J and Dctwcilcr S D Mora. Not. R. Astron. Soc. 211, 933 (1984) 
[86] Favata M, Hughes S A and Holz D E Astroph. J. Lett. 607, L5 (2004) 

[87] Tanaka T, Shibata M, Sasaki M, Tagoshi H and Nakamura T Prog. Theor. Phys. 90, 65 (1993) 

[88] Shibata M Phys. Rev. D 48, 663 (1993) 

[89] Papapetrou A Proc. R. Soc. London A209, 248 (1951); 

Dixon W G Proc. R. Soc. London A314, 499 (1970) 
[90] Tanaka T, Mino Y, Sasaki M and Shibata M Phys. Rev D 54, 3762 (1996); 

Suzuki S and Maoda K Phys. Rev. D 55, 4848 (1997); 

Suzuki S and Macda K Phys. Rev. D 58, 023005 (1998); 

Hartl M D Phys Rev D 67, 024005 (2003); 

Hartl M D Phys Rev D 67, 104023 (2003) 



